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ABSTRACT 

The probability density functions (pdfs) of molecular line centroid velocity fluctuations, and of line 
centroid velocity fluctuation differences at different spatial lags, are estimated for several nearby molecu- 
lar clouds with active internal star formation. The data consist of over 75,000 13 CO line profiles divided 
among twelve spatially and/or kinematically distinct regions. These regions range in size from less than 
to 1 more than 40pc and are all substantially supersonic, with centroid fluctuation Mach numbers rang- 
ing from about 1.5 to 7. The centroid pdfs are constructed using three different types of estimators. 
Although three regions (all in Mon R2) exhibit nearly Gaussian centroid pdfs, the other regions show 
strong evidence for non-Gaussian pdfs, often nearly exponential, with possible evidence for power law 
contributions in the far tails. Evidence for nearly exponential centroid pdfs in the neutral HI component 
of the ISM is also presented, based on older published data for optical absorption lines and HI emission 
and absorption lines. These strongly non-Gaussian pdfs disagree with the nearly-Gaussian behaviour 
found for incompressible turbulence (except possibly shear flow turbulence) and simulations of decay- 
ing mildly supersonic turbulence. Spatial images of the largest-magnitude centroid velocity differences 
for the star-forming regions appear less filamentary than predicted by decay simulations dominated by 
vortical interactions. No evidence for the scaling of difference pdf kurtosis with Reynolds number, as 
found in incompressible turbulence experiments and simulations, is found. We conclude that turbulence 
in both star-forming molecular clouds and diffuse HI regions involves physical processes which are not 
adequately captured by incompressible turbulence or by mildly supersonic decay simulations. The vari- 
ation with lag of the variance and kurtosis of the difference pdfs is presented as a constraint on future 
simulations, and we evaluate and discuss the implications of the large scale and Taylor scale Reynolds 
numbers for the regions studied here. 

Subject headings: ISM: clouds, ISM: molecules, ISM: kinematics and dynamics, stars: formation, 
turbulence 



1. INTRODUCTION 



tion of the vel ocity cor relatio n function and r elated 2-poiut 



a i,i i iji-ro-ii -u j iji. statistics (see Hobson 1992; Kitamura et al 
Although a great deal of effort has been devoted to quan = — i v c — j r c - r r- 



titativoly describing the complex column density spatial 



structure of star - form ing regions ( for recent approach es 3 ^e 



Langcr, Wilson fc Anderson 1993; Adams fe Wisemar 



Phillipj 1990, 19 91; | Cill fc Honrikoon l 1990 



1994; iWiseman fc Adams! 1994; E immcrmann & Stutzkil 
1993; IHoulahan fc Scalo| 1992; |Scalo| 1990; phappcll & 



Scalo] 1998; Stutzki ct al] 1998), comparatively little at- 



tention has been paid to characterising the radial velocity 
dimension of the data. Exceptions are studies o f possi- 
ble velocity dispersion -size s caling relations (see [ Falgar- 
onc, Pugct, fc Pcrault 1992, | Barranco fc Goodma n 1998, 



Goodman ot al. 1998, and references given there), ostima 



& Bally 1994 and references to earlier wor 



1993; |Miesch 
given there) 



an d searches for evidence of rotation (e.g. Goodman ct 
al.| 1993). Since one expects a signature of the dynami- 
cal and physical processes to appear in the velocity field, 
and because the velocity field is strongly coupled to, and 
may in a sense control, the density field and even th e star 
forma tion rate a nd the IMF (see models by, e.g. Fleck 
1983, |iik| 1995, |Scalo| 1998) it is important to develop 
additional diagnostics to investigate it. As a step toward 
a better understanding of mole cular cloud velo city struc- 
ture, Falga rone and co workers flFalgaronc 1989, Falgarone 



fc Phillips 1990, 1991, Falgarone et al. 1994; see below for 
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discussion) have explicitly tried to relate radial velocity in- 
formation to dynamical processes through the comparison 
of observed line profiles with frequency distributions, or 
probability distribution functions, found in experimental 
and simula tion studies of t urbulence. In an earlier short 
report we ( Miesch fc Scalc| 1995) extended that program 
to the frequency distribution of line centroid velocities. In 
the present paper we give a more detailed description and 
discussion of both the observational data and the centroid 
probability density functions (which we hereafter refer to 
as centroid pdfs) , and further extend the study to the pdfs 
of centroid velocity differences. Pdfs of line wid ths and 
line skewnesses are examined in a separate paper ( MieschJ 
Scalo |fc Ballyj 1999, hereafter Paper II). 

The primary goal of this work is to provide useful quan- 
titative observational constraints on ideas, models, and 
simulations of interstellar turbulence and its relation to 
star formation-constraints which must be accounted for 
by any theoretical approach that purports to provide a 
physical explanatory understanding of the phenomenon. 
A rec urring theme in the present paper is the degree to 



culled from early surveys of optical absorption lines and 
HI emission and absorption lines in order to show that 
good evidence already exists for roughly exponential cen- 
troid velocity pdfs in the lower-density atomic HI compo- 
nent of the ISM. The data employed in the present study 
of molecular star-forming regions are presented and dis- 
cussed in §3, and an overview of the dynamics, star for- 
mation activity, and physical environment in each region 
is provided in Appendix [A|. In §|] the statistical results are 
presented, including centroid velocity images and pdfs, as 
well as velocity difference images, pdfs and pdf moments. 
Several parametric and non-parametric pdf estimators are 
used and compared. For the velocity differences, the varia- 
tion of the pdfs with lag is emphasised and quantified and 
the possible relation of the non-Gaussian tails to filamen- 
tary structures is investigated. The centroid velocity and 
velocity difference maps are also used to derive effective 
Reynolds numbers and Taylor scales for the regions con- 
sidered. The results are summarised in This is the first 
in a series of papers on the statistical anal ysis of velocity 
fields in star- forming regions. In Paper II (Miesch, Scale 



which interstellar turbulence resembles incompressible tur- & Bally 1999), we use the same data sets presented here to 



bulence, a field with a very large experimental and theo- 
retical literature. However, we think it is important to 
recognise from the outset that even if interstellar turbu- 
lence does turn out to resemble incompressible turbulence 
in some respects, this would not imply that interstellar 
turbulence is "understood." Incompressible turbulence re- 
mains an essentially unsolved problem. As many authors 
have pointed out, there has been (arguably) little tangi- 
ble progress in the field of inc ompressible turbulence since 
Kolmogorov's seminal paper (Kolmogorov 1941, K41), and 
most of it consists of elucidating the ways in which K41 
was incorrect (e.g. the variation of the scaling exponents 
of different-order structure functions), and a proliferation 
of contrasting phenomenological theoretical approaches, 
mostly imported from other branches of physics (e.g. sta- 
tistical mech anical and field the oretical approaches; see 
She 1997 and L'vov & Procaccia 1997). If interstellar tur- 
bulence resembles incompressible turbulence, this result 
may help exclude some physical processes as being domi- 
nant, and may implicate the nonlinear advection operator 
in the momentum equation (which controls incompress- 
ible turbulence) as the dominant physical process. How- 
ever, the physics of this operator in the presence of a large 
number of degrees of freedom (i.e. large Reynolds number) 
is still obscure, and one must not confuse the existence of 
a huge literature with understanding. Furthermore, com- 
pressible turbulence does not entail any quadratic invari- 
ants (e.g. quantities conserved by the advection operator 
which are quadratic in the velocity, like kinetic energy), a 
property which is central to all models of incompressible 
turbulence. This should from the outset lead one to ex- 
pect that interstellar supersonic turbulence is different in 
fundamental ways from incompressible turbulence. From 
this point of view, the present work is aimed at uncovering 
observationally these differences. 

In §^ we present a discussion of the advantages and 
disadvantages of the use of centroid pdfs instead of line 
profiles, and a summary of previous work using this ap- 
proach (and the velocity difference pdf) for incompress- 
ible turbulence, molecular interstellar regions, and extra- 
galactic structure. We also present centroid velocity pdfs 



investigate spectral linewidth and line skewness variations, 
which provide yet another important diagnostic of molec- 
ular cloud velocity fields available from densely-sampled 
emission line observations. In the final paper of the series 
(Scalo, Chappell, & Miesch 1999, hereafter Paper III), we 
provide theoretical interpretations of the observed centroid 
velocity, velocity difference, and linewidth pdfs and a com- 
parison to analytic and numerical models. 

2. VELOCITY PDFS AND THEIR ESTIMATION 
2.1. General Background 

The one-point probability density function contains in- 
formation that is qualitatively different from correlation 
functions and related two-point statistics, which are mo- 
ments of some probability density. For example, the auto- 
correlation function or structure function, while containing 
spatial information, basically only involve the variance of 
the two-point probability density of velocities as a function 
of scale. The two-point probability density function itself 
could be computed, but is difficult to visualise, since, for a 
given lag scale, it is a two-dimensional surface, and its full 
representation would involve such surfaces at different spa- 
tial lags. More information about the velocity two-point 
pdf would require high-order moments, which cannot be 
constructed for ISM data because of noise. In contrast, 
the one-point velocity pdf, although it contains no spa- 
tial information, is easily displayed (it is basically just a 
one-dimensional histogram), giving access to the probabil- 
ity structure of the velocity field. Similarly, correlations 
between linewidth and other physical parameters only in- 
volve the variance of the one-point pdf, and thus averages 
over the information in the one-point pdf itself. For these 
reasons, the one-point velocity pdf offers a qualitatively 
different view of the velocity field. 

The relative independence of the velocity pdf and a sec- 
ond order moment of the two-point pdf, like the power 
sp ectrum or correlation function, c an be seen in the work 
of Dubinski, Narayan, & Phillips (1995), who generated 
simulated line profiles from random velocity fields with a 
prescribed power spectrum. Although they emphasised 
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the a bility of a Kolmogorov energy spectrum (k 5 / 3 ) to 



yield non-Gaussian line- profiles, their results show that 



were also found in hi gh latitude molecular clouds by BlitzJ 
Magnani, & Wandel (f 988). The presence of similar broad 



other forms of the energy spectrum would yield similar re- 
sults. In particular, a steeper k -2 spectrum, which might 
be expected for a field of discontinuities or shocks, also 
gave non-Gaussian profiles of similar form. Thus the 2- 
point velocity correlation function (a moment) is probably 
only weakly coupled to the 1-point velocity pdf. 

Recent work in several areas suggests that the one-point 
probability distribution function of dynamical variables 
like velocity is a useful tool that may be sensitive to dy- 
namical processes. These s tudies includ e large scale struc- 
ture of g alaxy velocities flBcrnardcai; 1994; Kofman ct 
al.Q|94; patelan fe Moscardini| 1994a,b), incompressible 
terrestrial turbulence (see below), distinguishi ng nonlin- 
ear chaotic processes from stochastic processes (|Wright fe| 
Schul{] 1993), and c haracterisation of s amples of musical 
volume fluctuations (Scalo & Chappell 1998). In particu- 
lar, studies of incompressible turbulence have shown that 
the higher moments (skewness, kurtosis,...) of the pdf can 
be used to constrain physical models for turbulent inter- 
mittency. Although non-zero skewness must exist at some 
level in order to provide energy transfer among different 
scales, the pdf of the velocity field itself in incompressible 
turbulence is in gen eral very n early Ga ussian, at least on 
large enough scales (Batchelor 1953 and Monin fc Yaglonj 
1971 review ea rly work; see mo re recent experim ents and 



simul ation s in Ansclmet et al. 1984, Figure 1; Kida & 



Mura ltami| 1989, F igure 6; |Jayesh fc Warhaft] 1991, Fig- 
ure 1; Chen et al. 1993, Figure 3.). An important ex- 
ception is the 3-dimensiona l incom pressible simulation of 
homogeneous shear flows by |Pumir (1996), who finds non- 
Gaussian, nearly exponential, velocity fluctuation pdfs for 
velocity compone nts perpendicular to the strcamwi se com- 



ponent. See also Lamballaise, Lesieur, & Matais (1997) 



for non-Gaussian velocity pdfs in channel flow close to the 
boundary. 

For incompressible turbulence non-Gaussian behaviour 
is well-established for velocity differences at small scales 
and velocity derivatives, and there is strong evidence from 
experiments and simulations for non-Gaussian behaviour 
in many other v ariabl es (see the papers referred to above 
and Chen ct al. 1989; Castaing, Gagnc fc Hopfinger 1990; 
Vincent fc Mcneguzzij 1991; [She ct"aL 1993). Often the pdf 



of the velocity difference or derivative field exhibits a near- 
exponential behaviour at smaller and smaller scales, and 
much work has gone into understanding this behaviour 
physically, especially in terms of the stretching proper- 
ties of the advection operator (see |5he| 1991 for a review). 
Part of the motivation of the present work is to investigate 
whether any of these properties occur in the more complex 
"turbulence" of interstellar clouds, and whether even the 
velocity fluctuation field itself presents measurable devia- 
tions from a Gaussian p df. 

Palgarone fc Phillips' (1990, 1991) have shown that line 
profiles constructed from high-sensitivity CO molecular 
line data (in several transitions) exhibit excess wing emis- 
sion, relative to a single Gaussian, over a very large range 
of scales, from 0.02 to 450 pc . For all these line profiles 
the width of the wings is about 3 times the width of the 
line core if both are fit by Gaussians, but the fractional 
intensity of the wing component (fraction of mass at high 
velocities) varies between about 0.03 and 0.8. Broad wings 



wings in regions whose scales are gravitating and non-self- 
gravitating, and in regions with and without internal mas- 
sive star formation, suggests that the behaviour is not due 
to stellar winds or collapse motions, and the variation in 
wing width in these regions seems to rule out a dilute 
warm gaseous component, as pointed out by Falgarone & 
Phillips. Since the line profile, in the optically thin case, is 
in effect a histogram of radial velocities, the broad wings 
have been viewed in the context of non-Gaussian pdfs, al- 
though there is some confusion concerning whether the 
line profile should be interpreted as the pdf of average line 
of sight velocities or of vel ocity differences; the latter in- 
terpretation is adopted by Falgarone & Phillips (1990) in 
comparisons with laboratory data. 

It is not clear that line profiles give a valid represen- 
tation of the velocity pdf. Every line profile samples a 
line-of-sight velocity field which in general contains a com- 
ponent whose characteristic scale is a significant fraction 
of the sample depth. The form of these systematic line-of- 
sight motions is unknown and may severely limit the corre- 
spondence between the line profile and velocity pdf. Such 
problems can largely be circumvented in analyses of simu- 
lations, where it is p ossible to to insure homogeneity o n the 
largest scales ( as in |Porter, Pouquct, fc Woodward 1994, 
as analysed by Falgarone et al, 1994), but homogeneity is 
probably not a good assumption in general for interstellar 
clouds. It is not difficult to show that the addition of a 
systematic component can significantly alter the estimate 
of the distribution of the velocity fluctuations, which is 
the function of interest. A cloud in non- uniform rotation 
about its center, for example, will yield non-Gaussian line 
profiles along lines of sight displaced from the projection 
of the rotation axis onto the plane of the sky (provided 
this projection is nonzero). In particular, these profiles 
will exhibit apparent excess wing emission due solely to 
the smearing arising from the variation of the line-of-sight 
component of the rotational velocity with depth in the 
cloud, which will thus distort the pdf of velocity fluctu- 
ations. In addition, radiative transfer effects can distort 
emission lines and cause the wings to become relatively 
more prominent if t he cloud is optically thick (although 



Falgarone & Phillips 1990 argue against this interpretation 



of the broad wings on the basis of their observed shapes) . 

An alternative procedure, which we adopt here, is 
to estimate the pdf of centroid line velocities (intensity 
weighted average velocity along the line of sight) sampled 
over a densely observed individual star formation region. 
While a "line profile" is a measure of the radial velocity 
(or velocity difference) pdf sampled along the line of sight, 
either in a single beam or averaged over many beams, the 
"centroid pdf" is the pdf of the mean velocity of line pro- 
files taken over a large spatial sample of positions in the 
plane of the sky. The two functions differ in the direction 
along which the sampling for the pdf is taken, and in the 
quantity sampled. 

The advantages of the centroid pdf approach include the 
much lower sensitivity required for each of the individual 
line profiles and the weaker dependence of the results on 
large scale systematic motions which, although still a con- 
cern, will tend to be mitigated by the line-of-sight averag- 
ing and by space filtering of the velocity maps (see §3.1). 
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For example, the centroid velocities of the rotating cloud 
discussed above will vary in an obvious way with position, 
and the effects of rotation can therefore be removed by 
applying an appropriate filter. Such a procedure is not 
possible with the individual profiles unless the rotation 
curve of the cloud is known. In addition, the presence of 
a warm "interclump" medium, or of "optical depth broad- 
ening", which would both contribute to the line profiles, 
will not much affect the pdf of centroid velocities, since 
the thermal component and the line saturation are sym- 
metric (although the centroid pdf, in the optically thick 
case, would only sample fluctuations on the leading edge 
of the cloud) . The problem with this approach is that the 
number of velocities (positions) which must be sampled in 
order to accurately estimate the tails of the pdf is very 
large, at least of order 1000. Furthermore, the relation- 
ship between the pdf of an average line-of-sight quantity 
(centroid velocity in this case) and the pdf of the radial 
velocity distribution in three dimensions is not clear. 

It is also possible to use the centroid velocities to con- 
struct the pdf of velocity differences for regions separated 
by a given distance, or, "lag" , and examine how this pdf 
depends on lag. These velocity difference pdfs are com- 
monly used in studies of incompressible turbulence (see 
references above) , and have been used to study models 



for th e cosmological evolut ion of galaxy velocities ( Peeblei 
1976; |Ueda, Itoh, & Satoj 1993; |Seto fc Yokoyamaj 1998J 



It is particularly intriguing that the galaxy velocity differ- 
ence pdfs exhibit exponential forms at small separations, 
very similar to the ISM pdfs reported here. 

A preliminary account of the velocity difference pdfs ob- 
tained for the star-formin g regions studied here was pre- 
sented by Miesch & Scalo (1994). Centroid velocity differ- 



ence pdfs 

region ( pis, Keene, et alj 1998; [Lis, Phillips, et al 



lave also been r eporte d for CO lines in the pOph 

" 1998) 



and f or HI lines in the Ursa Ma jor cirr us cloud ( Miville- 



Deschenes, Joncas, & Falgarone 1998). Lis et al\ (1996, 
1998) find that the velocity difference pdfs in simulations 
of mildly supersonic decaying turbulence exhibit strong 
non-Gaussian tails which are associated with filamentary 
structures and regions of large vorticity, and they report 
some evidence of such behaviour in molecular cloud ob- 
servations, at least in quiescent regions. Part of the work 
presented here is aimed at estimating the pdfs and spatial 
distribution of velocity differences for regions in which vig- 
orous high-mass star formation has taken place (Orion B 
and Mon R2 regions) although smaller regions containing 
only lower-mass YSOs are also represented (L1228, L1551, 
HH83). It should be recognised at the outset that the pdfs 
of the centroid velocities themselves for some of t hese re- 
gions already exhibit strongly non-Gaussian tails, ( Micsch| 
& Sca|k] 1995; see also §4.1 below), showing that the "tur- 
bulence" in star-forming molecular clouds is different from 
incompressible and mildly supersonic turbulence, which 
gener ally exhibit nea rly Gaussian velocity pdfs. Apart 
from Miesch & Scalo (1995), we know of only one other 



study in which centroid velocit y pdfs of molec ular line data 



have been presented, that of Padoan et al. (1997) who 
plotted histograms of centroid velocities for several areas 
in the Perseus cloud, another active star-forming region. 
Although large-scale motions were not filtered out and the 
pdf tails were not investigated in detail, these data too 
clearly exhibit non-Gaussian shapes. 



2.2. Previous Estimates for the Atomic ISM 

Work aimed at determining the pdf of interstellar gas 
motions dates back to the early 1950s. Several studies used 
optical absorption line velocities of "clouds" along the line 
of sight to OB stars and velocities of HI 21cm emission 
and absorption lines to estimate the peculiar velocity dis- 
tribution, after correction for solar motion and differential 
galactic rotation. These studies refer to fairly local gas, 
with distances less than about 500-1000 pc. However with 
the advent and subsequent prominence of molecular line 
observations, these studies were never repeated and were 
in effect forgotten. 

The results of these earlier studies are presented in Fig.l. 
In some cases the original histograms were only given in 
graphical form, and these plots were converted to digitised 
images using a scanner, and then measured on a computer 
terminal using standard image manipulation software. We 
estimate that the uncertainties in reproduction due to this 
procedure are at the 10 percent level. Since these stud- 
ies contain a small number of velocities compared to the 
present work, and because they are only meant to be il- 
lustrative of the results already in the literature long ago, 
the negative velocity portion of the pdf was reflected about 
v=0. 

Blaauw (1952) was apparently the first to estimate the 
form of the velocity pdf, using velocities of CaK lines from 



Adams' (1949) catalog. For 150 components along lines 



of sight to 120 stars with d<500 pc and 80 components 
toward 43 stars with d>500 pc, and excluding lines asso- 
ciated with the Orion region, Blaauw showed that in both 
cases there were too many components with high velocity 
to be fit by a single Gaussian, and that a single expo- 
nential could better fit the data. The probable existence 
of blends due to the relatively poor velocity resolution of 
the Adams survey was a problem, and Blaauw tried to 
correct his pdfs using an ingenious technique that is no- 
table for its use of humans as an analog random number 
generator. The resulting pdf for the 150 nearby compo- 
nents (Blaauw's Table 5a) is shown in Fig. 1 (filled cir- 
cles) in log-linear form. Notice that in this form a Gaus- 
sian would be a parabol a and an exponential would be a 
straight line. Takakubo (1967) also examined the K line 
velocities of Adams and concluded that an exponential pdf 
gave a better fit than a Gaussian. Huang! (1950) and Ka-| 
plan (1954) preferred 1/v fits to the optical line pdf , based 
again on Ad ams' catalogue. (See Fig. 1 in Kaplanj 1966). 



Siluk & Silk (1974) examined the pdf of the high-velocity 
(v > 20 km s^ 1 ) tail of the Adams optical line sam ple and 
fit the pdf with a power law of index around -3. [A4unch 



(1957) presented new observations of optical interstellar 
lines, at poorer resolution, toward 132 stars. He used the 
"doublet ratio method" (CaH/K or NaD2/Dl) to com- 
pare curves of growth with calculations based on Gaussian 
and exponential velocity pdfs. The Gaussian fit required 
a physically unrealistic increase of velocity dispersion with 
distance, while the exponential provided a fit with a single 
mean speed, confirming Blaauw's conclusion. 

Blending remained a problem with these interpr etations. 
Hobbs u sed inte rferometric scans of NaD (Hobbg 1969a, b) 
and KI (Hobbs 1974) in an attempt to fit individual pro- 
files, and showed that most of Adams' lines were probably 
multiple. Hobbs found that the majority of individual line 
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profiles could be fit by Gaussians, while a smaller fraction 
favoured an exponential, and many lines were not well-fit 
by ei ther distribution; however the velocity pdf was not 



exami ned. Falgarone fc Lequeux (1973) used Hobbs' NaD 



data to estimate the cloud peculiar velocity dispersion, but 
the number of components was too small for an estimate 
of the pdf. Later work on optical interstellar lines turned 
almost exclusively toward studies of interstellar depletion 
patterns, and the question of the velocity pdf was never 
re-examined. 



Fig. 1. — Histogram representation of the radial velocity probabil- 
ity density function (pdf) for distinct velocity components identified 
in optical absorption line surveys along the line of sight to OB stars 
and HI emission and absorption line surveys. The data were either 
taken from tables in the papers indicated, or measured from digi- 
tised scans of graphs presented in the papers. The negative-velocity 
portions of the pdfs have been reflected about v = 0. In this log- 
linear plot a Gaussian pdf would appear as a parabola while an 
exponential pdf would be a straight line. Additional studies not 
shown in this figure are discussed in the text. 

Meanwh ile, HI 21cm survey data was accumulating. 
Takakubc] (1967) performed Gaussian decomposition and 
subtraction of solar motion and galactic rotation, resulting 
in peculiar velocities for 544 HI emission lines whose veloc- 
ity pdf was presented for 3 intervals of linewidth and dif- 
ferent galactic latitude intervals. He judiciously refrained 
from any conclusions concerning the form of the pdf be- 
cause of the dependence of the derived velocity dispersions 
on galactic latitude. (This was perhaps the first paper to 
suggest consistency with the Kolmogorov velocity spec- 
trum for incompressible turbulence, on the basis of the 
power law scaling of velocity dispersion with the sine of the 
latitude corresponding to a scaling of the structure func- 
tion with path length.) The resulting histogram for all 544 
velocities is nevertheless shown i n Fig. 1 (reflected about 
v=0) as filled inverted triangles. Mast fc Goldstein] (1970) 
presented a similar analysis for 268 peculiar radial veloc- 
ities of HI 21 cm emission clouds that were well-resolved 
in velocity. Their Figure 4 clearly shows that an exponen- 
tial fits better than a Gaussian: a single Gaussian can fit 
the wings or the core but not both. Their observed pdf is 
shown in Fi g. 1 (open circles). 

Cr6visier (1978) studied the velocities of HI clouds iden- 
tified as components of Gaussian decomposition of absorp- 
tion profiles observed towards extragalactic sources in the 
Nancay survey. About 300 such clouds were selected on 
the basis of latitude, so that they might represent rela- 
tively local gas. The histogram of residual velocities, after 
corrections for solar motion and differential galactic ro- 
tation, is given in Crovisier's Fig. 1. The interpretation 
of the resulting histogram of residual velocities is prob- 
lematic because of the possibility that the intermediate 
velocity features are spurious. The histogram for all the 



components is displaye d in Fig. 1 (filled squares). Dickey, 



Salpeter, and Terzian (1978) presented the residual ve- 
locity distribution of clouds from the high-sensitivity and 
high-resolution Arecibo 21 cm emission/ absorption survey, 
but the number of velocities is too small to estimate the 



functional form of the pdf. Further discussion of the high- 
velocity tail of the HI distribution and implic ations for 
turbule nt energy r equirements can be fo und in Kulkarn: 
fc Fictj (1985) and [Lockman fc Gchman| (1991) 



All of the above studies included lines of sight that cover 
a significant fraction of the sky and a range in distance, 
rather than focusing on individual cloud complexes. An 
exception is the det ailed HI emission mapping of two re- 
gions by Verschuur (1974) using the NRAO 300 ft. an- 
tenna. Identifying "clouds" as distinct entities in spatial- 
radial velocity space, Verschuur presented histograms of 
velocities for about 150 clouds in his region A and about 
50 in region B. The region B histograms cannot be used to 
determine the intrinsic pdf of velocity fluctuations without 
careful filtering (which we have not attempted) because of 
ordered motion along the filamentary structures in the re- 
gion. The histogram for region A, reflected about v=0, is 
shown as open squares in Fig. 1. 

It should be emphasised that the rotation curve of the 
Galaxy and other large-scale trends have been subtracted 
out in all of the studies mentioned above. Furthermore, 
the results are likely not influenced by the Local Bubble, 
since the Verschuur sample was for a single region at a 
well-defined position and velocity, the optical-line clouds 
were detected toward OB stars which are on average about 
a kpc away, and the HI emission and absorption studies 
likewise sample much larger distances (the Local Bubble 
only extends about 100 pc from the sun). It is possible 
that the optical line studies are influenced by gas motions 
in local HII regions around the OB stars, so some caution 
should be taken in their interpretation. However, it is not 
clear that the clouds sampled in such studies lie in general 
at the same distance as the stars themselves (many studies 
have suggested otherwise), and the pdfs derived from the 
optical line data are similar to those derived from HI data, 
which are certainly not biased toward hot stars. 

Overall, it can be seen from Fig. 1 that evidence for 
non-Gaussian, and probably exponential, centroid velocity 
pdfs in the relatively low-density atomic HI component of 
the ISM has existed for some time. Many of these ear- 
lier works had already concluded that the velocity pdf is 
better fit by an exponential (or 1/v) than a Gaussian. It 
is surprising, and unfortunate, that all these early studies 
have been in effect "forgotten" in more current discussions 
of interstellar turbulence, since they suggest that at least 
some forms of interstellar turbulence differ significantly 
from incompressible turbulence, for which the velocity pdf 
is invariably nearly-Gaussian. The implication is that the 
physical processes involved are substantially different than 
for incompressible turbulence. This result for relatively 
low-density HI gas has recently received considerable sup- 
port from the d etailed aper ture synthesis HI study of the 



entire LMC by Kim et al. (1998). For scales above the 
resolution limit of 15pc, the peculiar (fluctuation) velocity 
pdf of the gas (for motions out of the plane of the LMC) 
is extremely well-fit by an exponential, and a Gaussian 
seems to be certainly excluded. 

We note that most of the HI studies reviewed in this 
section concern clouds which are probably not actively 
forming stars. So there is good evidence that, even in 
some "quiescent" regions, interstellar gas motions exhibit 
properties inconsistent with experiments and sim ulations 



of incompressible turbulence. On the other hand, Falgar 
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one et al. (1994) have found good agreement between the 
shapes (as quantified by moments of order 2 and 4) and 



50 kHz (0.14 km s" 1 ) for the L1551 data set. Typical 
signal-to-noise ratios are between 5 and 9. Further de- 



profiles in quiescent molecular clouds and synthcsiscd pro- 


be found in |B ally, Langcr, Stark, & Wilson (1987), Bally 


files fi 


om the "nearly incompressible" phase in simulations 


Stark, Wilson, & Henkel (1987), Bally, Stark, Wilson, & 


of mil 


ily supersonic, decaying turbulence by porter, Pou-| 


Langerl (1989), Pound, Bania, & Wilson (1990), Bally et 


quet, 


k Woodwardj (1994). Apparently, cither the physical 


al. (1991), Pound & Bally 1991 (1991), Bally et al.| (1995), 



processes occurring in at least some quiescent molecular 
clouds and HI clouds are substantially different, or the 
probability density function of centroid velocities contains 
different information than the second and fourth-order mo- 
ments of individual line profiles. The former posibility is 
likely, since the HI linewidths reported in the papers dis- 
cussed above are significantly supersonic in most cases. 

We also point out that shapes and shape variations 
in optically thin line profiles similar those reported by 



Falgarone et al. (1994) for mildly compressible turbu- 
lence simulations have also been claimed for c ompletely 
non-d ynamical incom pressible velocity fields (Dubinski, 
Nara)an, & Phillips 1995), highly sup ersonic simulations 
of cloud collisions (Keto & Lattanzio 1989), and highly 
supersonic stellar-driven self-gravitating turbulence fields 
by Vazquez-Semadeni et al. (Vazquez-Semadeni, private 
communication) . For this reason we are not convinced that 
the shapes and shape variations of individual line profiles 
alone can be used to discriminate between models, even 
for relatively quiescent clouds. Additional diagnostics of 
the velocity field are therefore needed. 

Thi s section has been a review of centroid velocity pdfs 
observed in the diffuse atomic phase of the interstellar 
medium. In the remainder of this paper, we extend these 
studies to the higher-density molecular gas in active star- 
forming regions and also consider the pdf of velocity dif- 
feren ces. Ultimately our goal is to use these statistics of 
the vc locity field to help understand the physical processes 



responsible for interstellar turbulence. 



3. OBSERVATIONS 



The observations used in this paper are the s ame as 
those us ed in the previous statisti cal analyses of M~iescb| 



fc Ballyl (1994, hereafter MB) and [Miesch fe Scale (1995, 
hereafter MS), and include large-scale mappings of the gi- 
ant molecular clouds in the northern part of the Orion 
region (Orion B) and in the Monoceros region, associated 
with the Mon R2 infrared cluster. Also included are ob- 
servations of the smaller-scale clouds which appear in the 
Lynds catalog as numbers 1228 and 1551 (Lynde 1962), 
as well as observations of the molecular gas surrounding 
the Herbig-Haro object HH83. In what follows, as in MS 
and in MB, we will refer to these five separate mappings 
as Orion B, Mon R2, L1228, L1551, and HH83. 

The Orion B, Mon R2, L1228, and L1551 observations 
are emission-line measurements of the 110 GHz, J — 1 — > 



and MB 

The observations of the HH83 molecular cloud were ob- 
tained with the IRAM-30m telescope on Pico Veleta, Spain 
in April and December 1989, and have a higher spatial res- 
olution (half power beamwidth = 13") and smaller field 
of view. The mapped region lies just west of the Orion A 
molecular cloud (L1641) and is associated with the Herbig- 
Haro object HH83. The measurements are of the 13 CO 
J = 2 — > 1 emission line at 220 GHz, with an effective 
frequency resolution of 100 kHz (0.13 km s _1 ). A typ- 
ical signal-to-noise ratio for these observations is about 
20. More information on the data acquisition, as well as 
the structure, environment, and physical conditions in the 



molecular cloud can be found in Bally, Castets, & Duvert 
(1994). 

We recognise that 13 CO is an imperfect tracer of the gas 
distribution and velocity field in dense molecular clouds. 
The primary problem is the likelihood that the 13 CO tran- 
sitions are optically thick along lines of sight with the 
largest total column density. Studies of the correlation of 
13 CO line strength with visual and inf rared extinction and 
with ra rer C O species such a s C 18 Q (Frcrking, Langer, fc 
Wilson| 1982; [Frcrking etH] 1989; [Lada et al] 1994) sug- 
gest that 13 CO saturates at visual extinctions greater than 
about 5 mag, although these papers only study a small 
number of individual clouds. Furthermore, the co mpari 



son of C 18 Q e mission with infrared extinction by Alves^ 
Ladgj (1998) suggests that for A v ^10 mag, C is O 



Lada, 



transition in 3 CO. They were obtained with the AT&T 
Bell Laboratories 7 m offset Cassegrain antenna located in 
Holmdel, New Jersey, between 1985 December and 1991 
June. This antenna has a very clean Gaussian beam with 
a FWHM of 100" between 98 and 115 GHz and the effec- 
tive channel bandwidth (resolution) is 100 kHz (0.27 km 
s _1 ) for the Orion B, Mon R2, and L1228 data sets and 

1 The centroid velocity in each line of sight is defined as T,TiVi/T,Ti, where T; and Vi are the brightness temperature and velocity corresponding 
to the i th spectrometer channel, and the summations span some selected velocity range. 



underestimates the total optical depth, possibly due to de- 
pletion of CO onto grain surfaces. These results imply that 
the present observations are biased in the sense that they 
do not penetrate the highest-column density substructure. 
However, the fractional area covered by column densities 
greater than A„ « 5 is small, about five or six percent in 
Orion B, one or two percent in Mon R2, and effectively 
zero in L1228, L1551, and HH83 (assuming the following 
conversion factors: I( 13 CO)^N( 13 CO) = 5.9xl0 14 cm" 2 
K" 1 km" 1 s, N( 13 CO)^N(H 2 ) = 7xl0 5 , and N(H 2 )^A„ 
= 1.06xl0~ 21 mag cm 2 ). So, we are probably probing 
most of the molecular gas in the surveyed regions. Fur- 
thermore, it would be extremely time-consuming to col- 
lect, say, C 18 spectra for the large number of lines of sight 
(more than 75,000) examined here, because of the required 
sensitivity. So, in terms of spatial coverage and statisti- 
cal characterisation of the velocity field, 13 CO is probably 
the best probe available. Nevertheless, the fact that this 
tracer probably does not probe the highest-column den- 
sity substructure should be kept in mind in interpreting 
the statistical results of §4 below. 

The observations are presented in Figures |and|. The 
coordinates corresponding to the origin (0,0) of each map 
are listed in Table [l] (see also MB, Table 1). The centroid 
velocity 1 along each line of sight is indicated by the color 
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Table 1 
Observed Clouds 



Region 


R.A. (1950) 


Decl. (1950) 


Distance a (pc) 


Size (pc) 


Mass a (M ) 


Mach Number b 


Orion B 
Mon R2 
L1228 
L1551 
HH83 


Q5 h 41 ™ 08 s . 1 
06 05 22.0 
20 58 00.0 

04 28 40.0 

05 31 06.3 


-01° 05' 00" 
-06 22 25 
+77 23 00 
+18 01 42 
-06 31 45 


415 [1] 
830 [2] 
300 [3] 
140 [4] 
480 [5] 


20x20^ 45x20 
45x40 
2.6x2.6 
2.1x1.5 
0.9x0.5 


8xl0 4 [6] 
9xl0 4 [6] 
200-900 [7] 
80 [8] 
8-15 [9] 


1.9-2.954.0-7.1 
2.6-3.9 
2.6 
1.7 
1.5 



References are as follows: [1] |Anth ony-Twarog| 198 2, [2] |Racine| 19 68, [3] |Greiiier ct al] 1989, [4] |Elias| 1978, [5] pcnzc] 
et al.[]l981, [6] [Maddalena et al.| 1986, [7] Pally et al] 1995, [8] |Snc!]| 1981, [9] |Bally, Castets, fc Duvert] l994 

b These Mach numbers are based on the centroid velocity fluctuation amplitudes listed in Table || (a c ), and an assumed 
kinetic temperature of 20K. 

c The first values listed correspond to regions la, lb, and lc, while the second values represent the remainder: regions 
2, 3, and 4 (see below). 



Table 2 
Regions Analysed 





Number 


Grid 


Interpolation 


Number 


LSR Velocity 


Region 


of 


Spacing 


Radius 


of 


Integration Limits 




Spectra 


(arcmin) 


(arcmin) 


Points 


(km s^ 1 ) 


Orion B (la) 


943 


0.5 


1.1 


822 


[9.35,11.75] 


(lb) 


1604 


0.5 


1.1 


2094 


[7.5,10.5] 


(lc) 


5817 


0.5 


1.1 


4801 


[8.5,12.5] 


(2) 


25510 


0.5 


1.1 


18420 


[5.0,14.5], [9.2,11.9], 












[7.0,12.1] 


(3) a 


23936 


0.5 


1.1 


10433 


[3.0,9.2], [1.0,7.5] 


(4) 


25311 


0.5 


1.1 


21876 


[7.5,12.7] 


Mon R2 (1) 


3260 


0.5 


1.1 


8630 


[10.6,13.4] 


(2) 


2813 


0.5 


1.1 


7592 


[7.7,12.5] 


(3) 


4570 


0.5 


1.1 


12208 


[11.2, 14.4] 


L1228 


1922 


0.5 


1.1 


3532 


[5.5,8.1] 


L1551 


2327 


0.5 


1.1 


4726 


[-9,-6] 


HH83 


1142 


0.05 


0.25 


6235 


[5.0,7.5] 



a As discussed in the text, region 3 in Orion B overlaps spatially with regions 2 and 4 and therefore includes most of 
the same spectra, but integrated over a different velocity range. 



image, with blue and red denoting movement respectively 
toward and away from the observer with respect to the 
mean LSR velocity of the cloud. The contour plot over- 
lays indicate the integrated intensity in the 13 CO transi- 
tion considered (see above). 

The number of independent spectra included in these 
observations totals more than 75,000 (see Table |^). The 
spectra for each region were extracted, reduced, and in- 
terpolated onto a regular spatial grid as described by MB. 
The grid spacing and interpolation radius used for each 



set of observations is listed in Table The velocity inte- 
gration ranges for each region, also listed in Table 0, were 
chosen to be large enough to span most of the emission, 
but still small enough to mitigate the influence of instru- 
mental noise on the centroid fluctuations. In the case of 
Orion B and Mon R2, some of the integration ranges were 
further restricted in order to isolate particular cloud com- 
ponents. See MB for further details. 

After interpolating the spectra onto regular spatial 
grids, threshold integrated intensity levels were introduced 
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Fig. 2. — Shown are the Orion B 3 CO J = 1 — » observations used in the analysis to follow. Contours represent integrated intensity 
and images represent centroid velocity. Red and blue in the centroid velocity images_correspond respectively to receding and approaching 
motions relative to the mean LSR velocity of the cloud (see the color table in Fig. ph. The data set has been subdivided into 6 distinct 
regions as described in the text. The orientation of the coordinate system with respect to the north and west directions on the plane of the 
sky is indicated. All coordinates represent the offset in arcminutes from the position of the origin, which is listed in Table hi and applies to all 
regions. The dashed blue line indicates the position of region 3, which overlaps regions 2 and 4. These, together with the aotted white lines, 
also serve to delimit the composite portions of regions 2 and 3 which were obtained using different integration limits as described in the text. 
The LSR velocity integration limits used for each region are indicated in brackets 



Fig. 3.— Similar to Figure |J but for the 13 CO J = 1 -> Mon R2, L1228, and L1551 observations and the 13 CO J = 2 -> 1 HH83 
observations. The color table used for the centroid velocity images (which applies also to Fig. H), is displayed on the right. The orientation 
of the coordinate system is indicated in each panel and all axes represent arcminute offsets with respect to the positions listed in Tabic hi 



in order to eliminate pixels with low signal-to-noise. Any 
remaining low-intensity emission located outside the main 
cloud boundaries and having dramatically different cen- 
troid velocity values (likely influenced by noise) was also 
removed. The number of pixels in each of the resulting 
centroid velocity maps in Figures | and | are listed in Ta- 
ble |^. Since the maps were oversampled, the number of 
points in the centroid and integrated intensity maps is in 
many cases larger than the associated number of indepen- 
dent spectra. 

Several distinct cloud components, or "regions" , distin- 
guished by their spatial position and their mean LSR ve- 
locity, emerge when the Orion B and Mon R2 data cubes 
are analysed closely. In the study presented here, as in MS 
and in MB, we have divided these data sets accordingly 
and have considered each of the components separately. 



Thus, in Figures 2| and 3a we have indicated six regions 
(la, lb, lc, 2, 3, and 4) in Orion B and 3 regions (1, 2, 
and 3) in Mon R2. In order to facilitate the isolation of 
specific cloud components in the Orion B and L1551 ob- 
servations, we have chosen coordinate systems which are 
tilted somewhat with respect to the right ascension - decli- 
nation standard. The orientation of each map is indicated 
by arrows in Figures |^ and ^|. 

All the coordinate axes in Figure || are with respect to 
the same origin in order to facilitate the placement of the 
different regions for the reader. Thus, region 1 (a, b, and 
c), lies directly above region 2 in the coordinate system 
chosen (see MB, Fig. le), or in other words, just north- 
east of the main Orion B cloud. It has been displaced and 
enlarged somewhat in Figure for clarity. Region 3 over- 
laps regions 2 and 4 in the manner indicated by the dashed 
blue lines. It appears to be a fairly well-defined separate 
cloud component which lies at a significantly lower LSR 
velocity (by ~ 4 km s _1 ) than the remaining emission 
and was isolated by choosing appropriate velocity inte- 
gration limits (see below and MB). It is likely associated 
with the doubly-peaked 12 CO line profiles in this region 



reported by Maddalena et al. (1986). In order to isolate 
the spatially overlapping regions in Orion B, we further 
subdivided regions 2 and 3 into portions, each integrated 
over a somewhat different velocity range, as indicated in 
Figure |^. These portions were then combined to produce 
composite maps for regions 2 and 3. 

The cloud components used for Mon R2 are shown in 
the three separate panels of Figure |^a. As in the case of 
Orion B, they were chosen on the basis of both spatial 
connectivity and mean LSR velocity using channel maps 
and spatial- velocity diagrams. The relative orientation of 
these regions can be discerned by noting the coordinate 
axes in each panel (which are all with respect to the same 
origin) and by referring to Figure Id of MB. The distinc- 
tion between the three different regions defined here is also 
apparent in the 12 CO channel maps presented by Xie & 



Goldsmith (1994; their Figures 1 and 2), who interpret the 
relative biueshift of the southeastern portion of the cloud 
(our regio n 1) as evidence for a large-scale, expanding shell 
(see §A.2 in Appendix A|). 

An overview of the physical conditions and environment 
in each of the molecular clouds in our study is given in 
Appendix |A|. Briefly, our sample covers a broad range of 
scales and conditions relevant to the velocity field, with 
energy input ranging from hot massive stars to low-mass 
protostellar outflow sources. What the regions have in 
common is that they all contain strong, active sources of 
momentum and energy to drive the "turbulence" whose 
statistics we wish to study. 

For reference, Table |l| lists estimates for the distance to 
each region, and for the total mass and plane-of-the-sky 
extent of the CO-emitting gas. Also listed are estimates for 
the Mach number of each region, based on the measured 
rms centroid velocity fluctuation (see Table || below) and 
on an assumed sound speed of 0.287 km s _1 (corresponding 
to a kinetic temperature of 20K). Note that the range is 
wide, from 1.5 to 7.1, and that all regions are substantially 
supersonic, suggesting that compressibility may play an 
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important role in their dynamics. 

4. STATISTICAL ANALYSIS 
4.1. Centroid Velocity Fluctuations 
4.1.1. Analysis 

It is apparent in the centroid velocity images shown in 
Figures g and || that many of the regions studied have 
large-scale velocity gradients across them, which at least 
in some cases, dominate any statistical quantifiers which 
may be computed. Since statistical approaches generally 
lose spatial information by considering only distributions 
or averages, their utility depends on some degree of unifor- 
mity in the data. Any statistical quantifier which is dom- 
inated by a small number of large-scale features in the 
data can give misleading results. In other words, large- 
scale gradients should be removed from the centroid data 
in order for our statistical results to be meaningful and 
comparable to theoretical models, numerical simulations, 
and laboratory experiments involving turbulent fluids. 

For these reasons, we have removed large-scale gradients 



by firsil applying a low-pass fillei, ui smoothing fmicliuii, 
to ea ch data set and then subtracting this smoothed map 
from the original to obtain the centroid velocity fluctua- 
tions. The filters were chosen to be as wide as possible, 
while still eliminating any prominent, extensive, antisym- 
metric "lobes" in the autocorrelation function of the fluc- 
tuation maps. Such features in the autocorrelation func- 
tion are the characteristic signature of large-scale gradi- 
ents in the data. The details of the filtering process are 
described by MB and will not be discussed further here. 

All of t he analy sis presented in this and the following 
sections (§ [l.l| -§ |4~3| ) is based on the centroid velocity fluc- 
tuations, which are shown as Grey-scale images in Figures 
U and ||. Although the pdfs of the centroid velocity fluc- 
tuations (discussed below) generally appear significantly 
different from the unfiltered data, they are not sensitive 
to the precise value of the effective cutoff wavelength or 
the shape of the filter. Fur thermore, since the analysis of 
velocity differences in §4.2 primarily samples small-scale 
structure, it is very insensitive to spatial filtering. 

Notice that although some of the highest-amplitude fluc- 



tuations, which will appear in the far tails of the corre- 
sponding pdfs, are found within the main bodies of the 
mapped regions, others occur near the cloud edges or in 
"outlier" regions-areas which could not be observed over 
large continuous spatial extents due to the sensitivity lim- 
its of the observations or the influence of overlapping cloud 
components. For example, this behaviour can be seen in 
the L1551 and Orion B, region lb maps in Figures |]c: and 
^p. Such high-amplitude edge or outlier points introduce a 
bias in the far tails of the estimated pdfs if they represent 
emission that is noise-dominated or kinematically distinct 
from the main body of each region. Although much of 
this emission was removed from the maps prior to analy- 
sis, it is very difficult to identify and remove all of it. The 
centroid velocity pdf estimates therefore still contain some 
of this bias. The influence of noise and overlapping cloud 
components is discussed further in Appendix [B| 

The pdfs, or probability density functions, correspond- 
ing to each of the centroid velocity fluctuation images of 
Figures | and | are shown in Figures ^ and [?]. Although 
the classical histogram is a familiar and robust est imator 
of the pdf, it is generally not the most accurate. Vio et 



al. (1994) have investigated several alternative pdf esti- 
mators in addition to the histogram and have compared 
their performance on some typical astronomical problems. 
They conclude that the histogram is decidedly inferior in 
the applications they considered and recommend the use 
of more sophisticated pdf estimators. Motivated by this 
result, we have im plemented two of the histogram alterna- 
tives described by Vio et al. (1994): the adaptive kernel 
estimator and the Johnson estimator. 

Kernel estimators approximate the pdf at each data 
value in terms of a sum over all data points, weighing 
each according to a specified window function, or kernel. 
They do not require binning of the data as in the classi- 
cal histogram estimator and are generally less sensitive to 
noise in the data, producing a smoother, continuous, and 
often more reliable result. In the adaptive kernel method, 
the effective width of the kernel increases systematically 
where the data density is low, thereby improving accuracy 
in the pdf tails. In our implementation, the form of the 
kernel and the manner in which its width is allowed to 



Fig. 4. — Shown are grey scale images of the centroid velocity fluctuations in Orion B, which were obtained by filtering the raw centroid 
velocity maps shown in Figure Has described in the text. Increasing brightness corresponds to a velocity directed increasingly away from the 
observer, as indicated by the intensity table in Fig. pi 



Fig. 5. — Similar to Figure but for the Mon_R2, L1228, L1551, and HH83 data sets. These maps were obtained by applying spatial filters 
to the centroid velocity maps shown in Figure H (see text). The sense of the grey scale is shown by the intensity bar, which also applies to 
Figure y. 
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vary are as described by Vio et al. (1994). The adaptive 
kernel estimation of each of the pdfs in Figures [| and [5] is 
shown as a dashed line. 

An alternative strategy for pdf estimation is based on 
defining some general class of analytic functions involv- 
ing several parameters which can be approximated using 
statistical moments or percentiles derived from the data 
set. The advantage of such an approach is that it does 
not require the subjective binning of data or the choice of 
a kernel function and window width. It can also be less 
sensitive to noise in the data than either histogram or ker- 
nel methods. The disadvantage of parametric methods is 
that, however general, they impose some functional form 
to the data which may not be present. Our implementa- 
tion of the Johnson system follows that described by |Vic| 
et al| (1994). Since the centroid velocity is in principle 
an unbounded variable, we consider only the so-called Su 
family of Johnson estimators, which is characterised by 
the four parameters 77, e, A, and 7, and which is defined 
by the following functional form: 



fj(v) = 



exp 



j 2 -K[{v-t)* + \*] 









7 + -q In | 




for — 00 < x < 00. The Johnson estimation for each of the 
pdfs in Figures I and I is shown as a solid line, and the 
best-fit parameters are listed in Table 

The centroid velocity pdfs presented by MS apply to 
the same data as the pdfs in Figures ^ and (?], but they 
were based solely on the histogram pdf estimator (with 
a slightly different binning than that used here). MS 
also constructed several composite pdfs, one of which was 
formed from the six subregions of Orion B, another from 
the three subregions of Mon R2, and a third composite pdf 
which was composed of all twelve data sets. For the sake 
of brevity, we have chosen to omit these composite pdfs 
from the present study, although it is worth noting that 
the Orion composite pdf in particular appears very nearly 
exponential. 



Most of the pdfs in Figures || and [?] exhibit significantly 
non-G aussian shapes, a result which will be discussed in 
$LL2| below. 



However, before proceeding, we note that 
a few pdf features probably originate from the superpo- 
sition of two or more distinct cloud components. For ex- 
ample, the morphology and kinematics of regions lb and 
lc in Orion B (see Figs. || and ||) suggest several cloud 
components, which may be the source of some of the bi- 
modality and fine structure exhibited by the associated 
pdfs (Figs, and ^|c). Although the spatially overlap- 
ping cloud components in Orion B regions 2 and 3 were 
separated out by choosing appropriate velocity integration 
limits, some contamination remains. Thus, the low veloc- 
ity (< —2 km s _1 ) pdf tail in region 2 (Fig. ||d) and the 
associated large kurtosis value of 11 (see Table || below), 
are probably influenced by residual emission from region 
3, and the abrupt high-velocity (~ 2 km s _1 ) cutoff in the 
region 3 pdf (Fig. He) is likely influenced by the imposed 
upper integration limit. Also, residual overlap in Mon R2 
may account for some of the excess in the high-velocity 
tail of region 1 and the low- velocity tail of region 2 (Figs. 
(1) |Tja and b). 

Another way to quantify the spread and shape of prob- 
ability density functions is by the sample moments, com- 
puted directly from the data set. We consider the first 



four central moments, defined as follows: 



1 

N 



%x, y v(x,y) ; 



(2) 



standard deviation = a = \j — S X)H [v(x, y) — /i] 2 ; (3) 



11 3 

skewness = — — [v(x, y) - fi] 
a iv 



kurtosis 



1 1 



^x, y [v(x,y) - /i] 4 



(4) 
(5) 



The summation over x,y spans the area of the cloud or 
in other words, the entire data set. The first two mo- 
ments quantify the location and spread of the pdf, while 
the third and fourth are nondimensional quantities which 
contain information about its shape. In particular, the 
skewness and kurtosis are measures of the symmetry and 



Fig. 6. — Centroid velocity pdfs for Orion B, corresponding to the images shown in Figure ]A. Crosses denote the histogram estimate with 
statistical (N 1 / 2 ) error bars and solid and dashed lines represent the Johnson and adaptive kernel estimators discussed in the text. The 
corresponding Johnson parameters are listed in Tabic tiL 



Fig. 7. — Similar to Figure but for the Mon R2, L1228, L1551, and HH83 data sets shown in Figure Q The parameters corresponding 
to the Johnson estimators (solid lines) are listed in Table M. 
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Table 3 

Johnson Parameters for the Centroid Velocity pdfs 



Region 


V 


e 


A 


7 


Orion B (la) 


1.55 


-1.63xl0~ 2 


0.192 


-0.114 


(lb) 


5.20 


-0.292 


1.01 


-0.292 


(lc) 


1.17 


0.218 


0.285 


0.504 


(2) 


2.03 


0.149 


0.534 


0.670 


(3) 


2.08 


-0.318 


1.13 


-0.604 


(4) 


1.32 


-7.75xl0~ 2 


0.476 


-0.157 


Mon R2 (1) 


5.37 


1.26 


1.46 


4.31 


(2) 


4.27 


0.121 


0.929 


0.609 


(3) 


2.80 


-6.70xl0" 2 


0.769 


-0.349 


L1228 


1.64 


5.00xl0- 2 


0.175 


0.322 


L1551 


1.11 


-5.65xl0- 2 


0.141 


-0.201 


HH83 


1.59 


5.29xl0~ 4 


0.158 


0.255 



"flatness" of the pdf respectively. A Gaussian distribution 
has a kurtosis of 3, and a value larger than 3 implies that 
the pdf in question has relatively more prominent tails - 
in other words, high-amplitude events are more numerous 
than would be expected for a Gaussian random variable. 
In the turbulence literature this behaviour is usually re- 
ferred to as "intermittency" (although the same term is 
often used in reference to the variation of the scaling ex- 
ponents of the moments of the velocity difference distribu- 
tion) . A kurtosis less than three implies the opposite-that 
pdf tails are less prominent relative to a Gaussian distri- 
bution. An exponential distribution {f(x) oc exp(— |a;|)) 
has a kurtosis equal to 6. Table | lists the second, third, 
and fourth moments for each of the pdfs shown in Figures 
|] and [t|. Note that the skewness and kurtosis values listed 
are central sample moments, and are therefore slightly dif- 
ferent from the values given by MS, which were computed 
using the histogram as a pdf estimator. 

Many theoretical and empirical studies involving pdfs 
in the context of incompressible turbulence are concerned 
with stretched exponential forms, w hereby the pdf tails 
fall o ff as f{x) oc exp {— a| x|^} (e.g. |Kailasnath, SrccnT 
vasan fc Stolovitzky|l992; |Lohse fc Grossmann| 1993; Fee] 
& Chjin|1997; prisch fc Somettej 1997). A Gaussian pdf is 
characterised by (3 = 2, and an exponential pdf by (3 = 1. 
Fractional (3 values can arise from random multiplicative 
processes and represent a straightforward and useful way 
of quantifying the departure of a pdf from Gaussian statis- 
tics. 

The value of the stretching exponent (3 determines how 
quickly the pdf decays, and as a result, how prominent 
the pdf tails are. For a properly normalised stretched ex- 
ponential pdf, there is a one-to-one correspondence be- 
tween the value of (3 and the pdf kurtosis. Therefore, if 
the observed velocity centroid pdfs can be approximately 
described by stretched exponential forms, then the sample 
kurtosis values listed in Table || imply an effective (3 for 
each data set. The results, listed in Table f|, yield stretch- 
ing exponents for Orion B and the smaller clouds (L1228, 



L1551, and HH83) which generally range between 1.0 and 
1.4. Exceptions include Orion B regions la ((3 — 2.1) and 
2 ((3 — 0.71). The Mon R2 regions yield higher (3 values, 
between 1.8 and 2. Curve fits to the tails of the Johnson 
pdf estimators yield similar, although slightly lower results 
for the stretching exponents. Adopting a Poisson weighing 
and combining the high and low velocity tails to improve 
statistics, such curve fits, excluding Mon R2, yield (3 val- 
ues ranging from 0.77 to 1.2, with the exception of Orion 
B, region lb ((3 — 1.6). The results for Mon R2 again 
indicate somewhat steeper tails, with 1.5 < (3 < 1.7. 

We have also considered the possibility that power law 
forms, f(x) ex x n , may fit at least some of the pdf tails 
better than stretched exponentials. Power law velocity 
distributions could arise from stellar winds, o utflows or 
superbubbles (e.g. |flk| 1995; |Oey fc Clarke) 1998) and 



would suggest physical processes substantially different 
than those occurring in isotropic, incompressible turbu- 
lence. When the pdf tails of Figures || and [?] are plotted 
on log- log axes rather than log- linear axes (Fig. ^|), many 
of them appear linear in portions, suggesting power laws. 
This is particularly the case for L1228, HH83 and Orion 
B, region 4, and possibly also in L1551 and Orion B, re- 
gions la, 2, and 3. Curve fits to the pdf estimators in 
these regions show some evidence for power law indices of 
k —4.5 ± 1, although the positive velocity tails in L1551 
and Orion B, region 3 seems significantly shallower (n w 2 
or 3). Power law forms are displayed in each of the frames 
of Fig. | for comparison. Although stretched exponentials 
in most (but not all; c.f. L1228) cases appear to provide 
somewhat better fits to the observed pdf tails than power 
laws, the data are generally not of high enough quality to 
make a reliable distinction. 

4.1.2. Discussion 

In most cases, there is good agreement between the three 
different pdf estimators (histogram, adaptive kernel, and 
Johnson) shown in Figures ^ and [?| Exceptions include 
the pdf excesses in the tails of Orion B regions lb and 2, 
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Table 4 

Moments and Stretching Exponents for the Centroid Velocity pdfs 



Region 


Standard Deviation 
(km s _1 ) 


Skewness 


Kurtosis 


P 


Orion B (la) 


0.15 


-0.018 


2.9 


2.1 


(lb) 


0.22 


-0.50 


4.8 


1.2 


(lc) 


0.36 


-0.58 


4.2 


1.4 


K 2 ) 


U.o4 


1 Q 
-1.6 


1 1 


n 7i 


(3) 


0.63 


0.073 


4.0 


1.4 


(4) 


0.46 


0.12 


4.9 


1.2 


Mon R2 (1) 


0.37 


-0.18 


3.0 


2.0 


(2) 


0.23 


-0.23 


3.2 


1.8 


(3) 


0.29 


0.14 


3.3 


1.8 


L1228 


0.14 


-0.17 


6.0 


1.0 


L1551 


0.18 


0.60 


4.6 


1.2 


HH83 


0.12 


-0.39 


4.2 


1.3 



Fig. 8. — The tails of the centroid velocity pdfs in Figures |fj and W are here shown on log-log axes in order to emphasise any power law 
components (oc v n ) which may be present. On such axes, powerJaw iorms would appear as straight lines. The abscissa in each panel runs 
from one to six times the pdf standard deviation, a (see Table W). Solid lines and asterisks denote the high velocity pdf tails and dashed 
lines and squares denote the low velocity pdf tails. All pdf tails are normalised to unity at a velocity of one a. Plot symbols (asterisks and 
squares) represent the histogram pdf estimator, while lines (solid and dashed) represent the adaptive kernel pdf estimator. The dotted line 
in the lower left corner of each panel illustrates a power law with a similar slope for comparison. 



and the fine structure in Orion B, region lc, all of which 
appear in the histogram and adaptive kernel estimators, 
but are not well-represented by the Johnson system. In 
general, however, the Johnson system provides a reason- 
able analytic approximation to the pdfs presented here, 
and the parameters listed in Table 0, together with the 
analytic form given by equation (|l|) can therefore be used 
by observers and modellers alike to roughly reproduce and 
compare with our results. 

The velocity pdfs for the Mon R2 subregions and the 
high-velocity portion of Orion B, region lb appear nearly 
parabolic on the log-linear axes shown, implying Gaussian 
forms, although the far tails for these regions still may 
suggest exponentials or power laws. However, most of the 
other pdfs appear distinctly non-Gaussian. Several of the 
pdf tails, particularly for regions 2, 3, and 4 in Orion B, 
as well as L1228, L1551, and HH83, instead suggest nearly 
exponential (straight lines on the log-linear axes of ^ and 
|^) , or possibly power-law (straight lines on the log- log axes 
of Figure ||) forms. These results are consistent with the 
derived stretching exponents, which lie near 2 for the three 
Mon R2 regions, and near 1.1 ±0.3 for most of the others 
(see Table and the associated discussion) . In no case do 
we find pditails which decay more rapidly than Gaussian 
(/3 <; 2 for all regions). 

It should be noted that the more Gaussian nature of the 



Mon R2 centroid velocity pdfs is not associated with par- 
ticularly high or particularly low Mach numbers relative 
to the nine other regions (Table nl) . However, relative to 
regions 2, 3, and 4 in Orion B, which have a comparable 
size and mass scale, the Mon R2 regions do possess some- 
what lower Mach numbers (Table and somewhat higher 
Reynolds numbers (Table 0). Also, differences in the dis- 
tribution and dynamical influence of young, massive stars 
in Orion B and Mon R2 may be significant (see Appendix 

0)- 

The pdf shapes can be described further using the kurto- 
sis values listed in Table [I| The three Mon R2 regions are 
all characterised by a kurtosis of about 3, again suggesting 
nearly Gaussian statistics. However, with the exception 
of Orion B, region la, the other pdfs have larger kurto- 
sis values, ranging from 4 to 6 (As explained above, the 
large kurtosis of 11 found for Orion B, region 2, is proba- 
bly anomalous, due to overlapping with region 3). These 
results are in strong contrast to those for incompressible 
turbulence, where th e ve locity pdf is nearly Gaussian (see 
references given in §2.1), and they imply a fundamental 
difference between the physics of interstellar supersonic, 
driven turbulence and incompressible turbulence. An in- 
terpretation of the difference based on the fact that kinetic 
energy is not a conserved quantity in supersonic turbulence 
is discussed in Paper III. 
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Our results for molecular gas powered by internal 
sources are basically in agreement with the velocity pdfs 
from optical line and HI emission and absorption compo- 
nents for cool atomic gas displayed in Fig. 1. Since these 
latter regions are generally not self-gravitating and contain 
no internal power sources, the similarity may imply that 
the non-Gaussian behaviour is either a result of nonlin- 
ear advection or that the latter regions possess turbulence 
that is powered by external sources. 

However, we note that at least 3 of the 11 subregions 
studied here do appear Gaussian, so our results should not 
be taken as an indication that all molecular clouds have 
exponential or even non-Gaussian velocity distributions. 
However since the 3 regions in question are all in the same 
GMC (Mon R2), we can say that the Gaussian distribu- 
tions are definitely in the minority. A similar statement 
can be made for the optical and HI line results discussed 
in sec. 2. 2. Perhaps the most interesting challenge, theo- 
retically, is to understand why the pdf is approximately 
exponential in most of the regions but apparently Gaus- 
sian in Mon R2. 

4.2. Centroid Velocity Differences 



As mentioned in §2.1, probability density functions of 



velocity derivatives or differences have proven to be an 
important and well-studied flow diagnostic in the context 
of incompressible turbulence, where they exhibit signifi- 
cantly non-Gaussian shapes. It is therefore of great in- 
terest to construct maps of centroid velocity differences in 
molecular clouds and investigate their statistical proper- 
ties. In this section we apply such an analysis to our data 
sets. 

The velocity differences we consider are defined as fol- 
lows; 

Aw(x, t) = i?(x) — u(x + t) , (6) 

where x and r are two-element vectors corresponding to 
position on the plane of the sky and v is the centroid 
velocity. The displacement t is referred to as the two- 
dimensional vector lag. The structure function is just the 
mean square of these velocity differences, averaged over 
the area of the cloud and normalised with respect to the 
centroid velocity variance (e.g. MB). In other words, the 
structure function at a particular lag is the normalised 
variance of the difference pdf corresponding to that lag. It 
therefore represents only part of the information available 
from the full pdf, which we consider below. 

Grey-scale images of the centroid velocity differences for 
each region are presented in Figures |^ and [l^, for a fixed 
vector lag, in pixels, of r = [1,1]. In other words, the 
images of Figures || and |l^ were obtained by displacing 
the maps of Figures ^ and ^ horizontally and vertically 
1 pixel and then computing the difference as expressed 
by equation (^). Of particular interest when considering 
the centroid velocity differences is the spatial distribution 
of large-amplitude velocity differences-those events that 
compose the pdf tails. Such events have been empha- 
sised in Figures || and |l0| by using a 3-level grey scale. 
As indicated by the intensity table in Figure [jlj, all pixels 
in which the magnitude of the velocity difference is more 
than twice the standard deviation of the map are shown as 
white. Magnitude values between one and two standard 
deviations are shown as light grey, while values less than 



one standard deviation are shown as dark grey. 

It should be noted that many of the extreme events 
found in Figures || and [l^ are likely spurious, influenced 
by low signal-to-noise ratios near the cloud edges and by 
residual emission from overlapping cloud components (see 
Appendix ^). This can be seen especially in Orion B, re- 
gion lb, L1551, and in the easternmost portions of Orion 
B, regions 2 and 4. However, the distribution of extreme 
events within the main volume of the clouds is more reli- 
able and generally reveals a very intermittent structure. In 
most cases, high-amplitude velocity differences appear to 
be well distributed throughout the clouds, although there 
is also some indication that some such events are corre- 
lated with the dense, star-forming cores of Orion B, region 
4 and Mon R2, region 1 (see also Figs. ||and||). At larger 
lags, the distribution of large-amplitude events becomes 
less intermittent, as demonstrated in Figure [ll] for Orion 
B, region 4. 

An important question is whether or not the largest ve- 
locity differences occur in correlated, filamentary struc- 
tures, as in the simulatio ns of mild ly supersonic, decaying 
turbulence discussed by Lis et al] (1996). We see little 
evidence for such structures in the images of Figures || 
and [l^, which appear to have a relatively more "spotty" 
spatial distribution. Although noise-induced fluctuations 
near low-intensity cloud edges probably do contribute, the 
relatively spotty appearance persists even in the cloud in- 
teriors and near those cloud edges (e.g. the intensity ridge 
in the western portion of Orion B, region 4) where the 
intensity is high and the influence of noise is small (see 
Appendix |b|) . 

The qualitatively different spatial distribution of these 
observed velocity differences is likely due to the fact that 
the regions discussed here are highly supersonic and are 
continually driven by stellar energy input. In this case the 
velocity field is expected to be dominated by the compres- 
sional (shock) modes, not the vortical modes that domi- 
nated the simulations discussed by Lis et al. It is impor- 
tant to recognizee that, except perhaps for scales smaller 
than about 0.1 pc, the turbulent velocity field of the ISM 
should be highly superson ic, and driven by shock interac- 
tions flKornrcich fc Scalo 1998). Our results for the ve- 
locity differences, along with the velocity pdfs themselves, 
as discussed in §4.1 above, therefore suggest a fundamen- 
tal difference between interstellar turbulence and incom- 
pressible turbulence, except perhaps in regions that have 
avoided energy input for a time long enough that they ap- 
proach incompressible conditions by means of turbulent 
decay. 

However Miville-Deschenes, Joncas, & Falgaron |(1998) 
report, based on their study of HI in the Ursa Major cirrus 
cloud, that the spatial positions corresponding to the non- 
Gaussian tails of the difference pdf are primarily located 
along three filaments in the cloud, although differences 
images like our Figures |^ and 10 are not shown. Such an 
image of velocity differences is presented in the CO (2- 



1) study of the pOph cloud by Lis, Keene, et al. (1998) 



Since the pOph core is a region with active internal star 
formation, it should be more directly comparable with the 
regions studied in the presen t paper. The ima ge of large 
velocity differences (Fig. 8 in Lis, Kccnc, et al. 1998) does 
indicate that some of the highest velocity differences are 
filamentary, although others are not. Thus the situation 
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Fig. 9. — Shown are grey-scale images of the centixud velocity differences for each subregion in Orion B, corresponding to a lag in pixels 
of t = [1,1]. As indicated by the legend in Figure hi the grey scale has only three shades, chosen to highlight the spatial distribution of 
large-amplitude difference values (white). 



Fig. 10. — Similar to Figure Bl but for the Mon R2, L1228, L1551, and HH83 data sets. Again, all maps correspond to a lag in pixels of 
■ = [1, 1]. The grey-scale intensity levels are based on the standard deviation of each map, as indicated by the legend. 



Fig. 11. — Shown are grey-scale images of the centroid velocity differences for Orion B, region 4 at the lags indicated (measured in pixelal 
The dashed boxes indicate the displacement as specified by the lag vector, t. The sense of the grey-scale is the same as in Figures H and llOl 



is not clear-cut, and the most we can claim is that, for the 
regions studied here, there is little evidence for filamen- 
tary clustering of regions with large velocity differences. 
Quan tification of this statement is beyond the scope of 



this p aper, since establishing the reality or absence of fil- 
amentary structure is a difficult and unsolved problem, as 
is well-known from studies of the large-scale distribution 
of galaxi es. 
As in § M4 

we can investigate the velocity differences 
further by considering their probability density functions. 
In order to facilitate visualisation and improve statistics, 
we combine data from all lag angles to obtain pdfs as a 
function of the lag magnitude, r = |t|. For each scalar lag, 
r, measured in pixels, there is an associated pdf which de- 
scribes the distribution of velocity difference values. The 
pdf variation with r is illustrated in Figures [l^ and [l3|. For 
each region, solid, dashed, and dotted lines correspond to 
scalar lags (t) of 1, 5, and 20 pixels. All pdfs have been 
compu ted using the adaptive kernel estimator described in 



4.1.1 



Some of the difference pdfs in Figs. |l2j and |13| appear 
to change from nearly-exponential at small lags to more 
nearly-Gaussian at large lags, a behaviour which is similar 
to incompressi ble turbu lence, the mildly supersonic decay 
simulati ons of [ Lis ct al.| (1 996), and t he observational re 



suits of Lis, Kcenc, ct ajj (1998) and Miville-Deschenes. 
Jonca s, fc Falgaron |(1998). However many of the regions 



show a persistence of exponential behaviour even at large 
lags (|t| = 20 pixels), although sometimes with a cutoff 
at very large velocity differences. This persistence of ex- 
ponential behaviour at large lags is seen in most of the 
observed regions, and again points to a difference between 
interstellar turbulence and incompressible turbulence. Al- 
ternatively, the absence of a clear transition to Gaussian 



behaviour may simply reflect the fact that, at the largest 
lags examined, the velocity field is still correlated while 
such larger-scale correlations do not occur in the incom- 



pressi ble regime o r in the regions examined by Lis, Keene 



et al. (1998) and N4ivillc-Dcschencs, Joncas, fc Falgaronc 
"(1998). 

The increase of the width of the velocity difference 
pdf with increasing lag is qualitatively similar to results 
found for incompressible turbulence. However the same 
behaviour may be expected for any stochastic field: at 
larger separations, the variables become less correlated, 
and hence a larger fraction of the differences have large 
values. Basically this result only shows that the structure 
function (variance of velocity difference) increases with in- 
creasing lag, a property shared by a large class of stochas- 
tic processes (e.g. Burgers turbulence, which is entirely 
compressible). In order to constrain models for interstel- 
lar turbulence, it is the functional form of this increase, as 
well as the behaviour of other moments, which is crucial. 

We therefore consider the variation with lag of the cen- 
tral sample moments of the difference pdfs, defined as in 
equations (|J)-(^|). The second moment, or variance, for 

as a function of the 
As mentioned above, 
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each region is presented in Figure 
scalar lag, r, expressed in parsecs. 
the variance of the difference pdfs is simply the struc- 
ture function, without the customary normalisation factor. 
The six subregions of Orion B are shown in frame (a) , and 
frame (b) shows the three subregions of Mon R2 together 
with the L1228, L1551, and HH83 results. The variation 
of the normalised fourth moment, or kurtosis [eqn. ([5 )1, 
with scalar lag for each region is displayed in Figure |l5. 

It is clear from Figures [lj and [l5] that the variation of 
the difference moments over some lag ranges can be ap- 
proximately described by a power-law dependence, which 
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Fig. 12. — Probability density functions of ccntroid velocity differences for the six subregions in Orion B. In each plot, solid lines correspond 



to a scalar lag, t, of 1 pixel, while dashed and dotted lines r ™ 



computed using the adaptive kernel estimator described in E 4.1.1 



•espond to scalar lags of 5 and 20 pixels respectively. All pdfs shown are 



Fig. 13. — Similar to Figure F^, but for the Mon R2, L1228, L1551, and HH83 data sets. Solid, dashed, and dotted lines again correspond 
to the velocity difference pdfs at scalar lags of 1, 5, and 20 pixels. 



appears linear on the log-log axes shown. A characteristic 
power-law index can be derived for each region by mea- 
suring the mean logarithmic slope of each curve over some 
chosen range in lag. However, the choice of the fitting 
range is not straightforward-most of the curves in Figures 
14 and [IB] exhibit steep slopes at small lags which flatten 
out as the lag increases. This "flattening out" is a common 
feature of all stochastic fields, which approach a constant 
variance and kurtosis at lags far exceeding the field's char- 
acteristic correlation length, or as in the present context, 
the effective width of the applied spatial filter. The fitting 
range should therefore not include large lags, where the 
slope is systematically small. On the other hand, sampling 
effects can influence the slope at the smallest lags (MB). 
We therefore chose the intermediate range of lags between 
6 and 12 pixels (note that this corresponds to a different 
physical scale in each region), and applied power-law fits 
of the form: 



variance oc r 



and 



kurtosis oc r Q4 . (7) 



The results yielded a 2 values between 0.64 and 1.05 (ex- 
cept for Orion B, region lb, in which a 2 = 0.33) and a 4 



values between -0.11 and -0.88, with typical values of 0.85 
and -0.5 respectively. These typical values are shown on 
the plots of Figures [li] and |l5f for comparison with the 
data. We emphasise that these best-fit power law ex- 
ponents depend significantly on the chosen fit range and 
should therefore be viewed only as a rough estimate of 
the logarithmic slope of the moment curves at intermedi- 
ate lags. With the exception of the L1228 variance («2 
= 0.66), the smaller-scale, low-mass star-forming regions 
L1228, L1551, and HH83 exhibit somewhat steeper log- 
arithmic slopes (1.01 < a 2 < 1.05, -0.88< a 4 < -0.64) 
than the GMC regions (0.33 < a 2 < 0.90, -0.52 < a 4 < 
-0.11), possibly reflecting the influence of molecular out- 
flows. The smallest slopes occur in the "cometary clouds" , 
regions la and lb of Orion B, although it should be noted 
that these regions have the poorest statistics (in terms of 
both signal-to-noise and number of spectra) . 

The variance results reported here are consistent with 
the more detailed structure function analysis of MB. The 
reader is referred to that paper for further work and for 
a comparison with structure function indices predicted 
by various analytic and numerical turbulence models and 



Fig. 14. — The variance of the difference pdfs is shown as a function of scalar lag, r, measured in parsecs. Frame (a) is comprised of the 
six subregions of Orion B and frame (b) includes the three subregions of Mon R2 and the remaining L1228, L1551, and HH83 data sets. The 
assumed distances to each region are as listed in Table hi. A power law of index 0.85 is shown in each plot for comparison. 



Fig. 15. — Similar to Figure |l^, but here the kurtosis of the difference pdfs is plotted as a function of scalar lag, t. As indicated by the 
legend in Figure [14, frame (a) corresponds to Orion B, and frame (b) corresponds to Mon R2, L1228, L1551, and HH83. A power law of 
index -0.5 is shown in each plot for comparison. 
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with empirical scaling relations from interstellar observa- 
tions (e.g. "Larson" scaling) and laboratory turbulence 
experiments (e.g. "Kolmogorov" scaling). 

The variation of kurtosis (flatness) with lag can be com- 
pared with the inco mpressible turbulence experiments of 



VanAtta fc Antonia (1980) and the multifractal models of 



Eggers fc Wang (1998), which exhibit a basically constant 



kurtosis out to some critical lag (which depends on the 
Reynolds number), beyond which the kurtosis decreases 
with increasing lag. However, except for a narrow tran- 
sition region near the critical lag, this decrease is much 
slower than that found in the present work. 

We believe that the ability of hydrodynamical simula- 
tions to account for the observed variations of second (Fig. 
14) and fourth (Fig. [l5|) moments of the velocity differ- 
ence pdfs with lag provides a decisive test of turbulence 
models, although it must be remembered that the obser- 
vations refer to regions that are supersonic and contain 
internal turbulent power sources. 

4.3. Reynolds Numbers 

The dynamical information available in the centroid ve- 
locity and velocity difference maps considered above can 
be used to assess the "degree" of turbulence in each re- 
gion, as quantified by the Reynolds number. The Reynolds 
number Re = vijv is defined only with respect to the 
scale £, which might be the size of the region, the cor- 
relation length (integral scale), or any other length that 
measures a useful characteristic spatial scale of the turbu- 
lence. One such scale which has proven especially useful in 
characterising incompressible turbulence experiments and 
simulations is the Taylor microscale, which measures the 
average spatial extent of velocity gradients. For incom- 
pressible turbulence, these gradients are entirely due to 
the vor ticity, so the Taylor microscale is defined as (e.g. 
Lesieur| 1990, p. 144) 



It = 



((V x u) 2 ) 1 / 2 



(8) 



where u rms is the root-mean-square fluid velocity. For 
compressible turbulence the velocity gradients are due to 
a combination of vortical (rotational, solenoidal) and com- 
pressible (irrotational, dilatational) modes. For highly su- 
personic turbulence the vorticity in eqn. (^) might be re- 
placed by the dilatation V • u. (In this case the simple 
scaling between Taylor scale Reynolds n umber a nd inte- 
gral scale Reynolds number [eq. VI-6-7 in Lesieur 1990, p. 
144] is lost because there are no relations between enstro- 
phy and viscous dissipation rate and the dissipation rate 
cannot be taken as u 3 /£ as for the dissipationless Kol- 
mogorov energy cascade.) 

The Taylor scale Reynolds number is generally smaller 
than the integral scale Reynolds number. For reference, 
simulations of incompressible turbulence are only capable 
of generating turbulence with Taylor scale Reynolds num- 
bers of at most 10 2 — 10 3 , while experiments reach Taylor 
scale R eynolds number of at most ~ 10 3 — 10 4 (s ee, for ex- 
ample , Kailasnath, Sreenivasan, & Stolovitzky 1992, 3hc 



et al. 1993, and references therein) 



For empirical studies such as the present one, it is ap- 
propriate to use the rms value of the measured velocity 
gradient in the denominator of eqn. (ph. This is clearly an 



important quantity, since it gives a measure of the char- 
acteristic scale over which turbulent interactions of any 
sort (e.g. shocks, vorticity stretching) occur. A problem 
is that estimation of centroid velocity gradients across the 
line of sight will be strongly amplified by uncertainties in 
the centroid velocities, while if we try to avoid this prob- 
lem by smearing the centroid velocity field with a filter, 
the resulting rms velocity gradient will decrease with in- 
creasing filter size. A convenient way to circumvent these 
problems is to take the appropriate mean velocity gradient 
as the square root of the variance of the velocity difference 
pdf taken at the smallest separation (one pixel), which we 
will call cr s , divided by the physical length corresponding 
to this separation, which we will call s. The reasoning is 
that the average taken when computing the variance of 
the difference pdf will greatly reduce the noise due to un- 
certainties in individual centroid velocities. We therefore 
define the Taylor scale as 



(9) 



where o~ c is the rms centroid velocity fluctuation ampli- 
tude. Taking the velocity occurring in the Reynolds num- 
ber as tr c , this gives a Taylor scale Reynolds number 



Rex 



cr 2 s 



(10) 



We take a c to be the dispersion of the mean spectral line 
profile for each region, given by the "parent dispersion" 
dp in Table 2 of MB. Alternatively, the "turbulent disper- 
sion", at in Table 2 of MB, could be used, but it includes 
uncertainties due to instrumental noise and it removes the 
influence of large-scale velocity variations under the as- 
sumption that they are uncorrelated with small-scale fluc- 
tuations. Since large scales in turbulent flows typically 
contain much of the kinetic energy and play an important 
role in the nonlinear energy transfer among modes, o~ p is 
the most appropriate dispersion measure to use in comput- 
ing an effective Reynolds number. However, the difference 
is relatively minor. 

Note that our definition of the Taylor scale (eqn. |^) 
and the associated Reynolds number (eqn. mi) depend on 
the scale at which they are computed-i.e. the spatial res- 
olution of the map. However, this is reasonable because 
it is the scale of the fluid motions we are sampling - we 
have no information about the dynamics on scales smaller 
than the resolution or larger than the emission region. The 
values may change if higher resolution maps become avail- 
able, but only if there is a large amount of power at high 
spatial frequencies. The range of values found for the Tay- 
lor scale, 3 to 10 times larger than the grid spacing, sug- 
gests that this may not be the case. Projection effects 
and line-of-sight averaging also influence the computation 
of the Taylor microscale, but these effects are very diffi- 
cult to quantify without some assumption about the three- 
dimensional nature of the velocity field. They could lead 
to an overestimate of the Taylor microscale and the cor- 
responding Reynolds number, although to our knowledge 
this has never been demonstrated. 

The molecular viscosity is estimated using a mean free 
path approximation, v = v t h/o~n, where Vth is the rms 
thermal velocity (kT//rni# ) 1 / 2 , n is the particle number 
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Table 5 

Reynolds Numbers and Taylor Scales 



Region 


L 


n 




s 




Re^/lO 7 


Rer/lO 6 


It 




(nc) 


()0 2 cm~ 3 ) 


(km s^ 1 ) 


(nc) 


(km s -1 ) 






(nc) 


Orion B (la) 


6.0 


2 


0.55 


0.12 


0.10 


1.2 


1.3 


0.66 


(lb) 


8.6 


2 


0.55 


0.12 


0.19 


1.7 


0.69 


0.35 




17 


2 


0.83 


0.12 


0.21 


5.1 


1.4 


0.47 


(2) 


30 


2 


1.15 


0.12 


0.21 


12 


2.7 


0.66 


(3) 


33 


2 


2.05 


0.12 


0.41 


21 


4.4 


0.60 


(1) 


40 


2 


1.31 


0.12 


0.25 


19 


3.0 


0.63 


Mon R2 (1) 


11 


6 


1.12 


0.24 


0.19 


50 


17 


1.4 


(2) 


32 


6 


0.74 


0.24 


0.12 


26 


12 


1.5 


(3) 


50 


6 


0.83 


0.24 


0.16 


15 


11 


1.2 


L1228 


2.6 


10 


0.74 


0.087 


0.064 


3.5 


13 


0.59 


L1551 


1.8 


10 


0.50 


0.041 


0.085 


1.6 


2.2 


0.24 


HH83 


0.67 


25 


0.42 


0.014 


0.042 


1.2 


2.6 


0.14 



density, and a is the collisional cross section. We adopt 
the cross section for H 2 — H 2 elastic collisions as about 
10~ 15 cm 2 . It is unclear to us whether the absence of a 
dipole moment for the H2 molecular should alter the elas- 
tic cross section compared to typical values for other atoms 
and molecules reported in the literature; we assume it does 
not. The particular densities are taken as the rough esti- 
mates given by MB: n(c m - 3 ) = 2500 for HH83, 1000 for 
L1228 and L1551, 600 for Mon R2, and 200 for the Orion 
B subregions. We adopted a characteristic temperature of 
20 K for all the regions. 

Then the Taylor scale Reynolds number can be ex- 
pressed as 



Re T = 1.8 x 10 6 n 2 cr 2 i5 Sp C /cr S! 5 



(11) 



where n2=n/100 cm , er c ,5 is the turbulent velocity dis- 
persion in km s , s pc is the linear scale corresponding to 
a one pixel lag at the adopted distance of the region, and 
17,5,5 is the standard deviation of the velocity difference pdf 
at this lag, in km s _1 . For comparison, the Reynolds num- 
ber corresponding to the scale L of the observation, taken 
as the geometrical mean of the semiminor and semimajor 
axis of the emitting region on the plane of the sky, is 



Re t = 1.8 X l0 b n 2 a C:5 L 



pc • 



(12) 



We computed these dimensionless numbers for each re- 



gion and the results are given in Table 



The region scale 



Reynolds number Rei, is, as expected, very large, vary 
ing between about lxlO 7 to 5x10 s . This rang e is very 



nearl y the s ame as the estimates of Re l given by Myers &j 
Khers jonsky (1995) for a sample of "diffuse," "dark," and 
"giant" clouds. Surprisingly, the Taylor scale Reynolds 
numbers, while generally smaller than Rer,, are still very 
large, between 7x 10 5 and 2x 10 7 . As mentioned above, the 
Taylor microscales, which measure the rms scale of veloc- 
ity gr adients, are, in units of the resolution s, a r Ja s —3 to 



10, suggesting that we are resolving the rms gradient scale. 
However, the very large values of the Taylor scale Reynolds 



number mean that, at the smallest resolvable scales, the 
advection term in the momentum equation, as measured 
by the velocity gradients, is still huge compared to viscous 
dissipation, implying that we have not yet resolved the dis- 
sipation scale, so velocity gradients must still be present 
on spatial scales smaller than the sampling grid of the ob- 
servations. Considering the numerical diffusion/artificial 
viscosity in existing numerical codes, simulations are far 
from approaching the resolution necessary to realistically 
represent interstellar turbulence. However, in the pres- 
ence of a magnetic field, the dissipation might be due to 
ion-neutral frictio n or magnetic reconne ction rather than 
ordinary viscosity ( Myers fe Khcrsonsky| 1995), making the 
situation less severe. 

A surprising result apparent in Table ^ is that the Taylor 
scale £t is approximately constant among most of the sub- 
regions of Orion B (£ T » 0.6 pc) and of Mon R2 (£ T » 1.4 
pc). The difference between Orion B and Mon R2 might 
be attributed to the differing spatial resolutions, since £t 
(eqn. ^ is proportional to the size of a resolution element, 
and Mon R2 is roughly twice as distant as Orion B. The 
same distance scaling is apparent for L1551, and qualita- 
tively for L1228 and HH83, where £ T is only 0.14 pc: HH83 
is at about the same distance as Orion B, but was observed 
at about eight times better resolution. These results may 
suggest that the Taylor microscale is approximately con- 
stant in all the regions observed when the reference spatial 
scale s (eqn. ^) is normalised to the same value. This is 
a particularly intriguing result, since it would imply that 
the characteristic rms scale of velocity gradients, measured 
relative to the rms global velocity field, is a characteristic 
length scale independent of differing physical conditions in 
the regions studied. However, we are unable to assign a 
value to this scale because it depends on the adopted scale 
over which the gradients are measured, taken here as the 
size of a resolution element. 



Experiments ( VanAtta fc Antonia 198 0; Tabeling ct al 
1996) and multifractal models ( |Bifcrak 1993; pggcrs fc 



Wang 1998) of incompressible turbulence indicate that 



the kurtosis (flatness) of the velocity difference pdf should 
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increa se with Re^ as about Ftj 15 , but for Taylor scale 
Reynolds numbers greater than about 700 the kurtosis de- 
creases with ReL- We find no evidence for such behaviour 
in the observations presented here, again suggesting that 
highly compressible turbulence differs significantly for in- 
compressible turbulence. However it should be noted that 
both ReL and Rer are proportional to adopted average 
densities, which are very uncertain, and so errors in the 
density estimates may mask any correlations. 

We have n ot attempted to estim ate magnetic Reynolds 
number (see Myers & Khersonsky 1995) for each region, 
because the magnetic viscosity depends on the magnetic 
field strength and the ionized fraction, quantities which are 
unknown for the present regions. The significance of these 
nu mbers for the propag ation of MHD waves is discussed 



Myers fc Khersonsky (1995) 



5. SUMMARY 

The probability distribution function (pdf) for fluctua- 
tions of molecular line centroid velocities, and line centroid 
velocity differences at different separations, have been es- 
timated for a number of local regions with active internal 
star formation. The data consist of a total of over 75,000 
13 CO line profiles covering five different molecular clouds. 
The internal stellar power sources include only low-mass 
protostellar winds in three regions (L1228, L1551, HH83) 
but extend to massive stars in the remaining two regions 
(Orion B and Mon R2). The GMCs Orion B and Mon R2 
were subdivided into six and three kinematically distinct 
subregions respectively in order to isolate particular cloud 
components. The total sample therefore numbers twelve 
distinct regions, each composed of about 1000 to 25000 
independent spectra. Centroid velocity fluctuation maps 
were constructed by interpolating these spectra onto reg- 
ular spatial grids and then applying spatial filters in order 
to remove large-scale gradients. These maps are displayed 
as images in Figs. |^ and pi. 

The pdf of centroid velocities for each region was es- 
timated using the classical histogram, a non-parametric 
adaptive kernel estimator, and the parametric Johnson 
est im ator. The methods generally agree well except in 



al. (1996). The pdfs are also markedly different from the 
nearly Gaussian velocity pdfs found in studies of incom- 
pressible turbulence. A theoretical interpretation of the 
centroid pdfs is presented in a separate paper (Paper III). 

We also constructed the pdfs of centroid velocity differ- 
ences for lines of sight separated by different scalar spa- 
tial lags. In agreement with other recent observational 
work on this function, and with studies of incompressible 
turbulence, we find nearly exponential difference pdfs at 
small lags which broaden with increasing lag. However, 
spatial images of the centroid differences show a "spotty" 
distribution for the largest velocity differences, with lit- 
tle evidence for the filamentary structures predicted by 
simulations of decaying marginally supersonic turbulence- 
structures which should be due to vorticity. 

We used the velocity difference pdfs to estimate the Tay- 
lor microscale, which is the rms scale of velocity gradi- 
ents. The Taylor microscale Reynolds number for all the 
regions is very large, ~10 5 to 10 6 , indicating that even at 
the smallest resolvable scales, the advection term in the 
momentum equation is still huge compared to viscous dis- 
sipation. 

Our observational results suggest fundamental differ- 
ences between turbulence in both star-forming regions 
and diffuse atomic interstellar clouds as compared with 
incompressible or mildly-supersonic decaying turbulence. 
These interstellar gas motions are characterised by very 
large Reynolds numbers, a high degree of compressibil- 
ity, continuous energy and momentum injection by inter- 
nal and external power sources, and various sources of 
anisotropy, including rotation, shear, magnetic fields and 
self-gravity. It might be necessary to incorporate some 
or all of these characteristics into numerical and theoreti- 
cal models before full agreement with observations can be 
achieved. The complexity of interstellar turbulence sug- 
gests that it may exhibit some degree of nonlinear self- 
organisation, producing coherent structures which could 
contribute non-Gaussian components to the velocity and 
velocity difference pdfs. Such behaviour has recently been 
demonstrated in three-dimensional numerical simulations 



the far tails of the pdfs. Although the Mon R2 regions ballaise, Lesieur 



of turbulent incompressible shear flow (Pumh 1996; Lam- 



exhibit nearly-Gaussian pdfs, except possibly in the far 
tails, all the other regions show strong excesses relative to 
a Gaussian, often suggesting nearly-exponential or power- 
law forms. Sample moments and stretched exponential fits 
to the pdfs were presented as quantitative measures of the 
departure from Gaussian statistics. These results confirm 
and extend the general conclusions reached by MS. 

Centroid velocity pdfs for diffuse interstellar HI regions, 
constructed from older published optical line and HI emis- 
sion and absorption line data, are also presented, and 
also show strong evidence for non-Gaussian, nearly ex- 
ponential, pdfs. The similarity with most of the molec- 
ular regions studied here suggests that either the expo- 
nential pdfs are not a product of stellar activity or that 
the internal velocity fields of the diffuse regions have been 
strongly affected by disturbances from external sources. 
These exponential centroid pdfs do not agree with the 
nearly-Gaussian pdfs found in numerical simulati ons of 
freely-decaying marginally supersonic turbulence by Lis ct 



ible convection (Brandenburg et al 



Matais 1997) an d turbulent compress- 
" 1996). On the other 



hand, it may be that some rather simple and fundamental 
property of the ISM is behind its statistical behaviour. For 
example, stellar winds, outflows, and super bubbles can all 
produce non-Gaus sian veloci ty distribution s with power- 
law pdf tails (e.g. |5k| 1995; |Oey fc Clarke] 1998). Alter- 
natively, as another example, the ISM is a system in which 
kinetic energy is not a global invariant conserved by the 
advection operator as it is for incompressible turbulence; 
in such a situation it is possible to obtain exponential ve- 
locity distributions, as arise i n simulations of s ystems of 
interacting wind-driven shells (Chappell & Scalo 1999). In 
Paper III, we will further address these and related issues 
regarding the theoretical implications and interpretation 
of the results reported in the present paper. 

This work was supported by NASA grant 5-3107 to JS, 
and by the National Research Council, through a postdoc- 
toral fellowship to MM. 



Miesch, Scalo & Bally 



19 



APPENDIX 



APPENDIX A: OVERVIEW OF REGIONS STUDIED 

All of the regions considered here are actively forming stars and there is evidence to suggest that the stars which have 
been produced have played a significant role on the dynamics of their parent clouds. In this appendix we provide a brief 
description of the environments and physical conditions in each of the regions studied in order to aid in the physical 
interpretation of the statistical analysis of §^ 

Orion B 

The Orion giant molecular cloud (GMC) complex is one of the nearest and most popular sites for studying massive 
star formation and its interaction with surrounding interstellar material. The GMC is comprised primarily of the A cloud 
(L1641), which lies behind the Trapezium stellar cluster and the well-known Orion nebula, and the B cloud (L1630), 
which is of a comparable size a nd mass and exte nds out to the northeast. The molecular gas in the region has been 



mapped e xtensi vely in 12 CO by |Maddalena et al 
Robinson! 1984 ; |Wood, Myers 
Bally * 



(1986), in the associated dust emission using IRAS (Beichman 1988; 



1991; Lada 



. ^ Daugherty 1994), an d in other isotopes of CO and CS (e.g. Bally et al 
Stark 1991; Kramer, Stutzki. fc Winncwisser 1996). The star formation activity and associated cloud dynamics 
have been reviewed by, e.g., Penzcl fc Stutzki (1989). 



The Orion GMC complex may have originated throug h the f ragmentation and local gravitational collapse of the ex- 
panding Gould's Belt supersh ell about 10 to 20 Myr ago Bally (1995) or from the collision of an infalling high-latitude 



cloud with the galactic plane (Franco et al. 1988). Massive star formation began by 12 Myr ago and produced the I Orion 
OB association. There is good evidence that these young, hot stars have had a profound influence on their parent clouds 
by accelerating, ionizing, ablatin g, co mpressing, and disrupting them throu gh supernov ae, rad iation, and stellar winds 



(Bally, Langer, Stark, fc Wilson 1987; Bally, Stark, Wilson, fc Langer 1989; Bally ct al. 1991; Gcnzcl fc Stutzki 1989) 



The large-scale velocity gradients and morphologies of the main Orion A and B clouds as well as the "wind-swept" or 
"cometary" appearance of a number of nearby smaller clouds, with "tails" directed away from the OB association (as in 
regions la and lb in the present study), provide a strong indication that the association has played an important role in 
their dynamics. 

The Orion B mapping considered here (Fig. |^) covers a number of regions which are still actively forming stars. The 
reflection nebulae NGC 2068 and 2071 are located near emission peaks in what we have defined as region 2, and the main 
peaks along the emission ridge in region 4 are near the well-known HII regions NGC 2023 and 2024 and the Horsehead 
Nebula B33. Other notable objects within the mapped region include a large number of refle ction nebulae, emb edded 



infrared s ource s, HII regions, Hcrbig -Haro objects, H2O masers, and molecular outflows (see Maddalena et al. 1986, 



Reipurth 1985, Lada, Bally, fc Stark 1992, Chandler & Carlstrom 1996 and references therein). In an infrared survey of 



the area, Lada, DePoy, et al. (1991) identified approximately 1000 sources, about half of which are probably embedded in 
the Orion B cloud. Most of these sources are clustered, with ~ 330 near NGC 2023 and 2024 in our region 4, and another 
300 near NGC 2068 and 2071 in our region 2. Active sta r formation is even occurring in at least some of the smaller 
clouds in region 1, as demonstrated by Reipurth fc Olbcrgj (1991), who report several Herbig-Haro jets and molecular 



outflows in L1617 (which we have labelled region lb). 

In summary, the Orion B observations studied here represent a very dynamic region where the molecular gas is sub- 
stantially influenced by supernovae, radiation, and winds from the massive stars in the nearby OB association and from 
embedded young stars which are continually forming within. 

Mon R2 

The Monoceros R2 molecular cloud is comparable in size and CO luminosity to Orion B. The l arge-scale distribu tion of 
12 CO follows well the area spanned by the Lynds dark clouds L1643, L1644, L1645, and L1646 flMaddalena et al.| 1986). 

The core of the Mon R2 cloud lies at the origin of the coordinate system used in Figure ||a, near the predominant CO 
emission peak, and is a well-studied region of ongoing high-mass star formation. It has been extensively mapped in CO 
isotopes, CS, HCN, H2CO, NH3, and HCO + , and is associated with a number of B stars and reflection nebulae, a large 
cluster of embe dded i nfrared sources, a com pact HII region, e nhanc ed X-ray emiss ion, a nd both H2O and OH maser s (see 
Bcckwith ct~al| 1976, fThronson et aT] 1980 , pontalban ct"at| 19 90, fTorrcllcs ctal] 1990, ponatas, Palmer, fc Novakl 1992, 



Giannakopoulou et al.| 1997, Tafalla et al 1997, Gregorio et al. 1998, and references therein). The embedded cluster has 
recently been studied in detail by Carpenter et al.| (1997), who estimate that its population numbers at least 475 stars. 



Also associated with this core region is one of the largest, most massive molec ular outflows known, exten ding for at least 
4 par secs a long its axis and involvin g almost 200 M & of interstellar material (Wolf, Lada, & Bally 1990; Meyers-Rice fcj 
Lada 1991; Xie, Goldsmith, & Patel 1993). The outflow is bipolar, centered near the middle of the embedded cluster, and 



dominates the dynamics of the molecular cloud within at least several arcminutes of the core (corresponding to th e origin 
in Fig ure ||a). A separate outflow, located about 75" from the main embedded cluster, was recently identified by Tafalla 



et al. (1997). 



A second prominent CO emission peak lies a bout 45' to the east (in region 2), and corresponds to the molecular core 



known as GGD12-15 (Little, Hcaton, & Dent 1990). Evidence indicates that this region too is actively forming stars 
and is also associated with a large bipolar molecular out flow, a compact HII region, several infrared sources, and an 



H2O maser (Rodriguez ct al. 1980; Olofsson & Koornnccf 1985; Little, Hcaton, fc Dcnt| 1990). Furthermore, there are a 
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number of IRAS point sources tracing the sharp emission ridges in regions 2 and 3 as well as additional infrared sources 
and reflectio n nebulae scattere d primarily throughout regions 1 and 2 and another HII region near the northeastern edge 
of region 1 ( Xie fc Goldsmith 1994). At l east 3 to 40 separat e infrared and X-ray sources have been identified outside 
the main Mon R2 core ( [Xie fc Goldsmith| 1994; |Gregorio et~aT| 1998). 

1 Loreil (1977) studied the large-scale kinematics of the region using CO, H2CO, and near-infrared observations and 
he interpreted velocity gradie nts across the clou d as a combination of rotation and collapse. However, the more recent, 
higher-resolution mappings by Xie fc Goldsmith (1994) have revealed a more complicated intensity and velocity structure, 
characterised by sharp emission ridges in the western and northeastern portions of the cloud and a southeastern portion 
with a lower LSR velocity and a more "wispy" appearance. This overall str ucture is confirme d by the observations 
presented here and is the basis for the region decomposition described in §0. Xie fc Goldsmith] (1994) interpreted the 
relative blueshift of region 1 as evidence for a large expanding shell, involving ~ 4 x 10 4 Af Q of material and moving toward 
us with a velocity of ~ 3-4 km s" 1 . Fairly recent, large-scale dynamical events are also suggested by the sharpness of the 
emission ridges in regions 2 and 3 and the evidence for associated star formation (see the preceding paragraph), which 
may indicate shock compression of the molecular gas. 

In summary, Mon R2 is a giant molecular cloud comparable in scale to Orion B. Although there is no analogous OB 
association disrupting the molecular gas to the same degree as in Orion, active star formation is indeed occurring within 
Mon R2 and there is also some evidence suggesting violent dynamical events on the scale of the cloud. 

LI 228 and L1551 

The Lynds dark clouds L1228 and L1551 are best known for the molecular outflows found within them. Both regions 
are nearby, well-studied sites of ongoing low-mass star formation. 

The L1228 cloud i s located in th e Cep heus flare, in a r ing of molecular gas which seems to be part of a 4xl0 4 year old 
supernova remnant (Grenier et al, 1989; Yonekura et al. 1997). At the center of the coordinate system in Figure Bp lies 
the young stellar object IRAS 20582+7724 an d associated with it , exte nding north east a nd southwest, is a large (18' x 
9'), w ell-collimated, bip olar molecular outflow (Haikala & Laureijs 1989; Bally et al. 1995; Anglada, Sepulveda, & Gomez 



1997; Tafalla & Myers 1997). The orie ntation and sense of the outflow axis is roughly the same as the large-scale velocity 
gradient apparent in Figure ||b. Also, Bally et al. (1995) have argued that there is at least one other YSO producing an 
outflow in the L1228 cloud core and have rep orted a numbe r of HH objects in the vicinity. Further evidence for active 
low-mass star formation has been provided by Ogura & Sato (1990), who identified 69 Ha emission stars in the area, nine 
of which are concentrated near the L1228 cloud core. 



L1551 occupies a small area along the edge of the Taurus- Auriga molecular cloud complex (e.g. Ungerechts & ThaddeuE 
1987) and is one of the most fertile centers of low-mass star formation known. The large (10' x 35'), collimated, bipolar 
molecular outflow in the southeastern portion of the cloud associated with the infrared source IRS 5 (which lies at the 
origin of the coordinate system in Fig. 0c) is one of the archetypical examples of the outflow p henomenon in young stellar 
objec ts (YSOs) and has been the subject o f much observational and theoretical re search (sec | Snell, Lorcn, fc Plambeck| 
1980, |Cabrit fc Bcrtout] 1986, |Uchida et al.| 1987, [Moriarty-Schieven fc Sneil| 1988, |Stocke et al.| 1988, [Bachiller, TafalTaJ 
fc Cernicharoj 1994~ Davis et al. 1995, and references therein). In the coordinate system chosen in Figure ||c, the outflow 
axis is oriented approximately horizontally, with blue-shifted gas toward the right of the origin and red-shifted gas toward 
the left, and it occupies much of the lower portion of the map. Also associated with IRS 5 and its accompanying outflow 
are a large reflect ion nebula known as HH102 and several other Herbig-Haro objects, including HH28 and HH29 (e.g. 
Graham fc Rubin] 1992). 

In addition to the IRS 5 vicinity, there are several other regions throughout the L1551 cloud which exhibit strong 
evidence for ongoing low-mass star formation. Especially several arcminutes to the north, near the very young T Tauri 
stars HL and XZ Tau and the YSO HH30, a large number of Herbig-Haro objects, mol ecular outflows, opt ical je ts, infrared 



source s, and excess Ha emission r egions have been identified ( Mundt, Ray, fc Buhrke 1988; Mundt et al. 1990; Graham & 



Heyer 1990; Pound fc Bally 1991 1991). Recent X-ray observations by Carkncr et al. (1996) and CCD and spectroscopic 
observations by Briceno et al.| (1998) have increased the total number of known T Tauri stars in the L1551 cloud to at 
least 26, along with one young B9 star. 

HH83 



The structure, dynamics, and nature of the HH83 molecular cloud have b een described in detail by Bally, Castets, & 
Duvert (1994). The innermost core has also been mapped in CS emission by Nakano et al. (1994). This cloud is smaller 
than the others in our study (see Table [l]), and is located at the extreme western edge of the Orion A cloud (L1641), within 
an area to the northwest of NGC1999 which has a high concentration of Herbig-Haro objects and several CO outflows, 
signifying active star formation in the region. 

An embedded star-the infrared source HH83 IR-lies in the southwestern part of the cloud, at the origin of the coordinate 
system shown in Figure ||d, and an associated optical jet extends o ut toward the northwest for at least 32", terminating 
in a conical bow-shock about 150" from the central IR source (e.g. Rcipurth 1989) . Also associated with th e embedded 
star and the jet is a low-velocity (« 5 km s" 1 ), poorly collimated, bipolar outflow ( [Bally, Castets, fc Duvert 1994). The 
outflow does appear in molecular emission lines, but there is good evidence that it is in a late stage of its evolution, having 
"blown out" of its parent cloud, and is currently depositing most of its energy and momentum into the atomic inter-cloud 
medium. Although the 13 CO-emitting gas does not trace the full extent of the outflow, small-scale velocity gradients near 
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the central IR source and evidence for evacuated cavities along the jet axis suggest that the outflow has had a substantial 
dynamical influence on the molecular gas in the vicinity. 

The axis of the outflow and jet, like several others in the region, is approximately aligned with the mean, large-scale 
magnetic field (directed northwest/southeast). The cloud is elongated perpendicular to this direction and, at least in the 
northern portion, exhibits a velocity gradient which may indicate gravitational collapse along the mean magnetic field 
lines on a timescale of ~ 3xl0 5 years. Also, the large-scale gradient across the cloud may indicate rotation about an axis 
aligned with the magnetic field and outflow. Furthermore, the cloud's proximity to the Orion OB association to the north 
suggests that it may have been accelerated southward relative to the surrounding gas by ablation, supernovae, and stellar 
winds. The associated compression may have triggered gravitational collapse in the southern end of the cloud, and may 
account for why star formation is currently occurring in the southern portion, but seems to be absent in the northern 
portion. 

APPENDIX B: SOURCES OF ERROR 



The influence of instrumental noise on the centroid velocity fluctuations was studied by |Miesch fc Baily) (1994), and 



we only summarise their results here. For each region studied in that (and in this) paper, they approximated the mean 
spectrum as a Gaussian with a peak, dispersion, and center determined respectively by the observed mean brightness 
temperature, lincwidth, and LSR velocity. The rms brightness temperature fluctuation due to instrumental noise was 
then estimated for each spectrum based on the effective integration times and detector temperatures, and a mean rms 
noise fluctuation was computed for each data set. A series of simulations was then performed for each region, based on the 
integration window used for the observations, the model Gaussian mean line profiles, and the mean noise fluctuations. In 
the simulations, white noise of the appropriate rms amplitude was added to the model profiles and the centroid velocity 
was computed. This operation was repeated (typically about 10000 times) with different random number seeds for the 
white noise, and the variance of the centroid velocity in these simulations was computed. The results indicate that the 
noise-induced fluctuations are very small, ty pically increasin g the standard deviation of the centroid velocity by only a 
few percent (compare a* to a c in Table 2 of Miesch fc Bally| 1994). 



Although these simulations do not directly address the influence of noise on the far tails of the centroid velocity 
and centroid velocity difference pdfs, they do suggest that this influence is minor, except for those spectra near cloud 
edges where the brightness temperature drops substantially below the mean brightness temperature. As discussed in §|[ 
threshold integrated intensity levels were used to eliminate spectra with very low signal-to-noise, but some questionable 
spectra near cloud edges still remain. For example, the extreme velocity and velocity difference values near the edges of 
Orion B region lb (Figs. | and §) and L1551 (Figs. | and 0) may be influenced by noise. 

In some areas, the influence of kinematically distinct, spatially overlapping emission may also bias the centroid pdfs, 
although the velocity ranges were carefully chosen to minimise this effect (see §H). This is particularly the case for the 
eastern portion of Orion B, region 4 (Figs. ^ and ^). 
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